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1. INTRODUCTION 

The principal objective of this report is to outline a multilevel control 
saieme for the Large Space Telescope (LST) . The concept and methodology of 
the scheme is based upon the deconposition-aggregation stability analysis of 
large-scale systems [1-3] , which was used to study structural properties of 
the control system for a spinning flexible spacecraft [4, 5] . 

The two -level analysis of the decomposition-aggregation method is ideally 
suitable for designing a multilevel feedback control [6-10] for dynamic systems 
conposed of interconnected subsystems. Local controllers on the subsystem 
level are used to stabilize (or optimize) the decoupled subsystems* On the 
second hierarchial level the global controllers are used to minimize the in- 
teractions among the subsystems, and make the control system meet whe required 
performance characteristics for the overall system. This multilevel strategy 
can solve complex control problems ”piece-by-piece" and make the conputer use 
attractive in cases when the direct approach is either not feasible Cexcessive 
conputer storage), or it is uneconomical Cexcessive conputer time). 

The detailed plan of the report is as follows: 

In Section 2, we will develop a nonlinear model for the LST which is based 
ipon the linear model described in [11] . The nonlinear representation will 
serve as a realistic model for evaluating the potentials of the multilevel 
schemes for control of the LST. 

In Section 3, we will outline the general multilevel stabilization algor- 
ithm [6-8], Both local and global controllers are involved. The local con- 
trollers are used to stabilize each decoupled subsystem by any of the classi- 
cal techniques such as pole- shifting, root-locus, parameter plane, etc. The 
role of the global controllers is to minimize the effect of interactions among 
the subsystems* Finally, the aggregate system is constructed on the higher 
hierarchial level to conclude stability of the overall conposite system. ‘We 
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will consider a class of dynamic systems |12] idiidi can always be stabilized 
by the proposed scheme using local controllers only. Since the LSI model de- 
veloped in Section 2 is in that class, we will be able to effectively design 
the feedback control which stabilizes the LST, 

In Section 4, we will present a multilevel optimization scheme for con- 
trol of large-scale systems [9, 10]. The local controllers are lased to opti- 
mize the decoi:^led subsystems with respect to quadratic cost. The global con- 
trollers are applied to reduce the subsystem interactions, or entirely decouple 
the subsystems as is the case of the LST. 'While this control scheme results 
in a suboptimal performance when the effective interactions are present, it 
produces an optimal control when the total decoupling takes place. Thus, the 
design procedure can effectively be used for constructing an optimal control 
system for the LST, 

Both the stabilization and the optimization multilevel schemes are en- 
tirely conputerized. The description of the programs is provided in the Ap- 
pendix, 

This report is written under the supervision and with the participation 
of the Principal Investigator, D. D. Siljak. Investigator S. K. Sundareshan 
developed the model of LST in Section 2, and Sections 4 and A. 2 on multi- 
level optimization. Investigator M. B. Vuixevic developed the multilevel sta- 
bilization scheme presented in Sections 3 and A.l, 
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2. DEVELOH«r OF A MODEL FOR H-IE LST 

The Large Space Telescope (XST) is modeled as a rigid body witli tliree 
orthogonally mounted reaction wheel actuators and is considered to be subject 
to gravitational and magnetic disturbance torques. Unlike in the earlier ana- 
lyses [11] , nonlinear coupling phenomena are not ignored and a con^lete tliree- 
axes model for the spacecraft is obtained as a nonlinear interconnected sys- 
tem. Tlte interconnections represent the coupling between the motions along 
the individual axes. Hence this model mil be a more accurate description of 
the LST, which however, is necessary due to the precision pointing requirements 
demanded of the control system. 

Hie spacecraft’s equation of moticii can be mitten down from the Euler 
equations [11] , as 

I • w + u X I » « + y {w X w tr I. + 2w. X I- • u + I* " «J* 

* .'I IL X X X XX 

1=1 

+ (ij. X I. ■> no.) = M (2.1] 

111 ^ 

and 

I. • <Jj. + u) X I. » w + CO X u. tr I. + 2w. X I. * to 
11 1 1111 

+ I. * £ 0 . + to- X I, • to* = M- , i = 1, 2, 3 (2.2] 

iiiiii' 

where 

I is the inertia tensor of the LST given by 
I = 

I , I , I denoting the components along the three axes 
*c z 

constituting an inertial reference frame Irf! 

, i = 1, 2, 3, are the inertia tensors of the tliree reaction 
wheels that are mounted orthogonally and parallel to the 


10 0 
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axes constituting the standard body-fixed reference frame 
®rf and hence can be expressed as 


ii = 


^Ix ^ 


"ly 

0 




0 0 



0 ; 


I 


2z 


I 


3 




0 0 



M is the angular velocity vector of the LST relative to the 
frame ; 

, i = 1, 2, 3, are the angular velocity vectors of the re- 
action wheels relative to the frame i^£ ; 

M is the total external torque acting on the LST; and 

, i = 1, 2, 3, are the internal torques on the reaction 
wheels , 


The angular velocity to can be expressed in terms of the rates of angu- 
lar deviations along the three axes of as 



to = 


0 



C2.3) 


where is the roll angle, 0 is the pitch angle and iJj is the yaw angle. 
Similarly the angular velocities to^ of the reaction wheels can be expressed 




ill torms o£ tlie cojnponents as 



C2.4) 


Equations C2.1D and (2.2) can now be simplified into the following four 
sets of scalar equations: 

(i) Equations governing tlie motion of the LST body: 

Ijc? + sfl'Cij- y + I32V3 5 - ^1x^1 ' “x ] 

V * ■ Hz’s ^ “ ’’V 

+ iecy I3.) + l2yV2 * ■ ^ * ^Si.^3 = ^ J 


C2.S) 


(ii) Equations for the Reaction hheel mounted parallel to x-axis: 

= ''ax 

V ^ V®* * Vi* “ V ■ 

V 'V = hz , 

(iii) Equations for tlie Reaction fcel mounted parallel to y-axis 

^2z^* ^^2y“ hy^2^ ^ ^*2x 

^2y® hy^2 " ^*2y 

hz^ ^ - ^2yV = 

(iv) Equations for the Reaction hlieel mounted parallel to z-axis; 


(2.7) 


1 : 1. - J . - I i - .J- 

6 

For further simplification, we will assume that the reaction wheels are 
small so that , l^y « 1 ^ , I 32 « degree of 

freedom ■'nly. With tliese, equations C2.53’‘C2.8) can be sinplifxed into, 


y * eici,- ly) + ■ 


V * * hyh ° 

* 

ij.? + l 0 cy y + 132^3 = ^ 

C 2 . 9 ) 


( 2 . 10 ) 

hy^ * hy'’Z “ “27 

( 2 . 11 ) 

hz^ * ° “ 3 z 

( 2 . 12 ) 


Substitution of (2.10) -(2.12) into (2.9) will result in the folloidng 
three equations describing the motions along tlie individual axes and their 
interconnections ; 


y . a,- y = (m^- ' 

lyS + Cy = 5«y- M2y} ■ 

* cy y ei = (M^- M3P ^ 


C2.13) 


It is now necessary to evaluate the various torques. Since the internal 
torques on the reaction wheels are small, it may be assumed that these are 
proportional to the control signals actuating the wheels. Hence, 


^ 2 y " "^^2 ^ 
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Wh^jre 

in (2,14) merely indicate the directions of these torques). 

The external torques acting on the body of the LST are mainly environment- 
al disturbance forces and are con^osed of gravity- gradient, magnetic, aerody- 
namic and solar pressure torques. The latter two mil be negligibly small com- 
pared to the others and mil usually be accounted for in control system designs 
by considering them as equivalent zero-mean stationary white noise processes. 
The gravity-gradient and magnetic torques can be represented as purely deter- 
ministic signals involving a constant term and a sinusoidal function of time 
with twice orbital rate . Hence , following the analysis in 111] the external 
torques can be obtained as. 




are the drive motor constants (.the negative signs 



^^11 ^12 
^^21 *^22 
^^31 *^32 


cos(ait + x) + s^} 

cos( t + x) + s„> I 
^ 7 

C0S( t + X) 'S’ s,} I 
i) Z 




(2.15) 


Wliere Yj^j » i = 1, 2, 3 , are constants that can be determined [31] from the 
inertia conponents I^, 1^, , the magnitude of the LST dipole moment and the 

earth's magnetic field intensity; and s^, i = 1, 2, 3 , are white-noise pro- 
cesses characterizing the aerodynamic and solar pressure torques. 

Substitution of (2.14) and (2.15) in (2.13) and further sinplification re- 
sults in the following system of equations: 


= 3^u^ + 

0 + - $ 2^2 ■** ^ 
ift + ag'jB « PjU^ + 


(2.16) 


tdiere a. 


. V 


0,n = 




3 I, 





1 




“ and M , M , M are the external disturbance torques given by (2.15). 
X y z 

It is now sintple to obtain a state-space representation o£ the LST by 
choosing tlie state-vector 


E‘1', 5>, 0, 0, ifl , 


C2.173 


which results in the time-invariant model, 


where 


X == Ax + Bu + h-Cx) + FM 


C2.183 


0 1 0 0 0 0 
0 0 0 0 0 0 
0 0 0 1 0 0 
0 0 0 0 0 0 
0 0 0 0 0 1 
0 0 0 0 0 0 


03^ 0 


^2 


0 3. 


0 0 
0 0 


The diagonal structure of the matrices A , B and F permits us to par- 
tition the state -vector as, 


X = [xj^,X2,X3] 


C2.19) 


where 


^ ^1 ^ 


=^21 _ ® 


=^31 

n ' "^3" 

0 Xto 


g 


With, this, (2.18) can be described as a set of interconnected subsystems, 
■vdiere 


0 

L_ 

11 


0 

1 


1 

0 

1 

0 

0 

, b. = 
’ 1 

fi. 

II 

1 

1 



0 


0 


0 

= 




, hjCx] ^ 



^*^1^2^32 _ 


j;“2^32^2_ 




and dj = being the external disturbances. 

It may be observed that when hj^C^xO =0 ? i = 1, 2, 3 , (2.20) represents 
three decoupled subsystems •vdiich describe the motions of the spacecraft along 
the three axes. However, h^(^) ^ot zero and constitute the interconnec- 
tions among the subsystems, thus maldng an analysis based on the smaller -dimen- 
sional decoupled subsystems alone inaccurate. 

The system represented by [2.20), is driven by the disturbance forces d^ 
in addition to the control signals . However, these external disturbances 
can be con^letely cancelled by constructing a disturbance accommodating control- 
ler as described in [3, 1]. This involves the determination of a suitable dif- 
ferential equation model for the disturbances and with the augmentation of the 
disturbance variables with the state variables of the system, designing a feed- 
back controller that counteracts the disturbance forces by feeding back the es- 
timated disturbance variables. Although this analysis is conducted for a single- 
axis model of the LST (only for the pitch motion control] in [1] , a straight- 
forward extension that iises three separate disturbance accommodating controllers 
can he obtained for the three -ax5,s model presently considered. Due to the above 
reason, we will ignore the disturbance terms from our model and conduct all fur- 


with d^ = M^, d 2 = M^. 
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ther analysis on the system, 

x. = A.x* + b.u. + h. Cx] » i 2, 3 , (2.21) 

X XI IX ' 

obtained from (2.20) with the substitution = 0 . 
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3. STABILIZATION 

IVheii a complex dynamic systCTi i? given, as a number of locally controSled 
interconnected subsystems, it can be stabilized by a multilevel control sclieme 
[ 6 - 8 ] based upon the decon^osition-aggregation stability analysis 1 1 - 3 ] . 
In the scheme, the dimensionality problem is resolved by carrying out all op- 
erations on the subsystem level. Both local and global controllers can be in- 
volved. The local controllers are introduced to stabilize each decot^jled sub- 
system by any of the classical techniques such as the pole-shifting by state 
feedback, root-locus, parameter plane method, etc. The global controllers 
minimize the effect of interaction'’* among the subsystems. Finally, the aggre- 
gate system is constructed on the liigher hierarchical level to conclude stabi- 
lity of the overall composite system. 

It is important to note that the proposed stabilization produces large- 
systems which are connectively stable [1 - 3 ] . That is, stability is in- 

variant under structural perturbations whereby subsystems are disconnected and 
connected again in various ways during tlie operation of the system. Further- 
more, the stabilized systems liave wide tolerance to nonlinearities in the in- 
teractions between the subsystems. 

After we outline the multilevel control scheme for stabilization of large- 
scale systejns, we mi] consider a class of dynamic systems which can he always 
stabilized by the scheme using local controllers only. Since the LST model 
developed in the preceeding section falls in that class, we mil be able to 
effectively design the feedback control which stabilizes the LST. 

3.2, Ccntrot 

Let us consider a linear dynamic system 

X “ Ax + Bu , C3.1!) 
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where xCt) G is the state of the system, uft) 6 R^ is the input to the 
system, and A and B are constant n x n and n x s matrices. We assume 
that the system is brought into the input-decentralized form 

s 


^i 


= A-x. + 
1 1 


.1 


A. -x- 
3-J J 


■fb.Ui 


i = 1, 2, . 


C3.2) 


n. 

1 


where G R is the state of the i-th subsystem, and 6 R is 

the corresponding local control, so that 


r“=R-^x R‘x ... xR 


A 


n. 


n 


(3.3) 


and each pair (A^, b^) is controllable. 

In (3,2), the matrices A^, A^^ , and the vectors b^|^ have appropriate 
dimensions. As shorn in reference [ 6 ], any linear dynamic system (3,1) 
can be represented by its input-decentralized form (3.2). 

To stabilize the system (3.2), we apply the decentralized feedback con- 
trol 

u^(t) = u^(t) + u|(t) (3.4) 


where u^(t) is chosen as a local control law 

"i = •’'K > 

n- 

with a constant vector k^ G R , and u|(t) is chosen as a global control 
law 



S rp 

I hT.x. , 
j=l ^ ^ 
j?^i 


(3.6) 


n. 

where k- . G R ^ are constant vectors. 

By substituting the control (3,4) into (3,2), we get the closed- loop sys- 
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tern as 


X. = 




- 1, 2, 


( 3.73 


Since each, pair (A^, bj) is controllable, a single choice of can 

T 

be always made [13] to place the eigenvalues of A^- at any desired 

* ■ « « ■ * ■ 
locations ± jaij, ... , ± 3“p, . ... j “»J;, (cf^ > 0 ; q = 1, 2, 

... , , and 1 ^ p ^ . Then, each. uncoiq)led sybsystem 


— (A^** t i~l>2, ... ,s 


C3.8) 


is stabilized with a degree of e^qionential stability 


IT. = mm CT 
^ a ^ 


C3.9] 


To provide a Liapunov function [ 5 - 8 ] with the exact estimate of 
for each decoupled subsystem, we apply to (3.8] the linear nonsingular trans- 
formation 


^i “ W » 

to get the system (3.83 9-s 

% = , 

-1 T 

where A. = T. (A.- b.kOT. has the quasidiagonal form 

•L »i> 3 b 3b IL 


C3.103 


C3.ll) 


■ 

i 

i 


i 


• 



“1 


-a 

P 

CO 

P 

i ^ t 

< 


-4 

, . 0 ■ f 

i 

“CO 

i 

»"^p4*l* *■“ ’ "V [ 


L 1 

i_ 


L p 

PJ 



(3.12) 
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n. 


For ths system (3.11), ive choose the Liapunov function v^: R ^ ^ Rj. ? 

C3.13) 




where 


and 


T~ 

K *n 


A-K. + H.A. = "G. , 
11 11 1 ' 


(3.14) 


^i ^ ^®i ^^1» **' ’ °p+l’ '** » ®i^i * 

(3.15) 

In (3.15), ^ ^ is an arbitrary constant and is the n^^ x n- 

identity matrix* 

The aggregate comparison system involving tlie vector Liapunov function 
v: -s- R® , 


J- 1 1 


^ V2» ... , v^) , 


(3.16) 


is obtained for the transformed system (3.7), 
” A.x. + y (A.." b-kT.)x- , 

•i 11 ij X ly j ^ 


X. 


jj4i 


1 ' 1, 2, ... , s 


(3.17) 


where A- - = tT V. .T. , b. == T^^b. , RT. = kT.T. , and using the Liapunov 

functions ‘V‘^(x^) defined in (3.13). Using the aggregation method presented 
in [1 ” 5 ] , we construct the conparison system 

v<W, ( 3 vl 8 ) 


where the constant s x s matrix W = (w. -) has the elements defined as 


W* • = -5.. IT- + (1”6..}£-- . 
Ij 13 1 '' 1T^^13 ' 


13 ■'"13 


(3.19) 


where 6.. is the Kronecker symboi, ir. is defined in (3.9), and 

Ij T X 
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C3.20) 


■where is the inaximum, eigenvalue o£ the indicated jnatrix. 

As known [1-3] , global asyn^totic stability of the system (3.17} and, 
therefore, original system (3.2), is implied by the Sevastyanov-Kotelyanski con- 
ditions [14 ] , which for W = (w.j) defined by (3.19) and (3 ..20} have the 
following form 



“'^1 

^12 •*' ^Ik 



j 

^21 

““^2 ^2k 

>0, k = 1, 2, , s . 

(3.21) 


% 

** 

hz *’* "’^k 


i 

i 


To satisify conditions (3.21), we choose -the vectors k.. in (3.20} so as 

Ij 

to minimize the nonnegative numbers which reflect the strength of inter- 

connections among the subsystems in (3.17). Such choice is pro'vided by 


k.. = [(bfb.) •*- b-A. 
ij 1 r 1 ID-* 


(3.22) 


where (h^b^) b^ is •the Moore-Penrose generalized inverse of b^ [15] . 

■fe M A 

The choice of . in (3,22) produces the optimal aggregate ma-trix W in 

XJ 

A 

the sense that W* (that is, W - W 0} is valid for all , That 
is equivalent to saying [16] that 


conditions (3. 21} are necessary and sufficient for 0 , that is, for 




stability of W , the choice k^^ = k^^ is justified. 


To conclude stability of the overall system (3,17) with the optimal choice 


-* 


k^ j = k^ j , which is 




j = 1, 2, . . ' 


(3.23) 




1 
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we apply the determinantal inequalities C3.21) to the optijnal aggregate matrix 
= CWy) defioied by (3.19] and = 5?^. = AjJ 

We arrive at the following: 


T'heoi?Bm 3,1, The Hnecw contool system (3.2] is stabitised by the linear con- 
trol taws 

s 


"i = ■’‘k ■ \ ' i = 1. 2, 


, S 


(3.24] 


T -*T -1 

where k*. = k. -T. , if the corresponding s x s aggregate matrios 

J 

«* = ♦ Cl-5i3)lyl 


(3.2S3 


satisfies conditions (3.21]. 


Successful application of the above theorem depends on appropriate choice 
of the eigenvalues for the decoupled subsystems (3.8]. Once the subsystem 
eigenvalues are prescrioed, the control law (3,24] and, thus, the gain vectors 
k^, k^^. in (3.24], are cou^juted uniquely using the proposed algorithm. There- 
fore, if for Gonputed gains k. , k.. , the conditions (3,21] are not met, a re- 
assignment of the subsystems eigenvalues is required. The search for aiy ap- 
propriate set of subsystems eigenvalues can be aided by the interactive con^u- 
ter program described in Section A.l. The efficiency of the computer program 
relies on the low order of the subsystems and the sinplicity in testing the 
Sevas tyanov-Kotelyanskii conditions (3.21], Furthermore, the conputerized 
procedure provides a considerable freedom to the designer to apply his under- 
standing of the system and the familiarity with the method to ccme with a 
successful design. 
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3*3 An Ittustrat'UDa Example 

Let us consider a system C3.1) described by the equation 


1 

11.50 

86.50 

0.45 

0 

-4.09 

0.18 

1 

8.36 

0 

1.25 

14.75 

0.18 

0 

-9.82 


The eigenvalues o£ the matrix A 


4 

22.50^ 


■"1 

0 

8,91 

-0,82 


1 

0 

0,36 

3.27 

X + 

0 

0 

5 

2o7S 


0 

1 

0.18 

-6.36_ 


_0 

1_ 


corresponding to (3.26} are 


C3.263 


and the system is unstable. 

To stabilize the system (3.26), ive start with its input-decentralized re- 
presentation (3,2) given as 


^1 = 


' 1 11,50 86.50' 

0.45 Q -4.09 


Xi-f 


' 4 22.50' 

8,91 -0.82 


''2 ^ 


1 

I 



Lo.lS 

1 

8,36j 

Lo.36 

3.27J 

1..C 


" 5 

2,75'" 


0 

1.25 

14,75“ 


"I*' 

= 



x„ + 







_0.18 

-6.36__ 


_0.18 

0 

-9.82_ 

1 

_1_ 




u, 


2 * 


(3.28) 


and transform each subsystem into its comparison form [13] to get 


^‘1 



0 

1 

0“ 


~ 3.20 

1.98“ 


rn 


0 

0 

1 

4 * 

-14.72 

0.49 

4 * 

0 


^8.86 

8.50 

9.36 


_-7.92 

36,01__ 

_1_ 





0 1 " 


“ 1,69 

1.26 

0.08“ 


“o'" 

ft 


X, + 




Xt + 



__32.32 -1.36 _ 

ti 

__-7.52 

-S.23 

0.49 J 

X 

^1_ 


«2 C3.293 
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The trajisfoxmtion into the conoparison form is of no conceptual signifi- 
cance, and is pexfomed on, the subsystem level for two practical reasons. 
First, it is convenient for subsystem stabilization by pole -assignment apply- 
ing the state feedback and, secondly, the diagonal form C3.17] with no complex 
roots, can be obtained from the con^anion form C3.29} using the Vandermonde 
matrix in (3.10) where x^ is replaced by x? . 

Now, by using the local feedback law C3.5) and vectors 

k^ (1791.14, 458.50, 46.36] 

kg == (33.82, 1.14] , (3.30] 

we allocate the eigenvalues of tlie uncoupled subsystems (3.27] from 


to 


1 " > 

^ 5.04 


,g -1.39 , 

= -6.41 


>2 - 10.12 , 


(3.31] 

1 = -10 . 

X^ = -1 


Z ” ’ 

4 = -1-5 


1 = -15 . 


(3.32] 


After the local stabilization, the interconnected subsystems have the quasi' 
diagonal form 


• 

“10 

0 

0" 


"-23.52 

-43. 7S" 


“0.1 


0 

-12 

0 


40.23 

68.97 

^2 

-0.17 


_ 0 

0 

-1S_ 

^-15.49 

-24.96 _ 


„0.07_ 


« ' 

“rl 

0 


179.95 

247.53 

367. .4.0"“ 


"2" 

•w 

^2" 

_ 0 

-l.S_ 

^2 

-182,58 

-249,03 

-365.95 J 

Xi + 

-2_ 


(3.33] 
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wMch. is not identical to C3»17), For the moment, we did not mahe use o£ tlie 
global control ^ C3.33}. In order to demonstrate the effect of tlie 

global controllers, we set 1^2 ~ ^21 ~ ^ * 

From C5.9) and (3.32), we have = 10 t "^2 ~ C3.203 and (3.33] 

we compute |^2 “ 98.51, - 676.68. The aggregate matrix W in (3,18] is 

obtained as 


’ -10 98.51 

676,68 -1 


(3,34] 


which does not satisfy the conditions (3.21], Therefore, we cannot conclude 
stability of tlie overall system. 

Let us use now tlie global control specified by (3.22], 


= (-238.95, -415,34] 

= (90.63, 124.14, 183.33] (3.35] 


whidi yields the subsystems (3.33] as 


• 

lio 

0 

0” 


"0.37 

-2.25“ 

A* 

0 

-12 

0 

^1 

0.40 

-0.25 


- 0 

0 

-15 _ 


_ 0.44 

2.73_ 


1 

0 

X, + 

"1.31 

-0.75 

0.72“ 

0 

-1.5_ 


^1.31 

-0.75 

0.72_ 


(3.36] 


and the aggregate matrix 


-10 3.55 

2.37 -1 


(3.37] 


•which satisfies the conditions (3.21]. Hierefore, by ■dieorem 1 the system 
(3.28). is stabilized by the contToI law (3,24] determined by the gains (3.30] 
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and C3.3S). The eigenvalues of tlie overall closed-loop system 


-10 

0 

0 ! 0.37 

-2.25“ 



0 

-12 

0 i 0,40 

-0,25 



0 

0 

-15 i 0.44 
1 

2.73 

X 

(3.38) 







-1.31 

-0.75 

0.72 I -1 

0 



-1.31 

-0.75 

0.72 1 0 

-1.5„ 




corresponding to (3.36), are 

2 “ - jO.16, = -10.27, « -11.99, « -15,17, (3.39) 


■Mhidi have negative real parts. 

It is also interesting to note that an upper estimate of the degree tr 
of es^onential stability of tlie syston (3.1) is provided by the aggregate 

« A 

matrix W since, in general tt < min ir- . In other words, the degree of 

” i ^ 

exponential stability of the overall system tt stabilized by the proposed 
method, is smaller than the degree of exponential stability of each decoupled 
subsystem. 


3.3 Local Stabilization 


In this section, we consider a class of linear input- decentralized large' 
scale systems which can always be stabilized by only local feedback control 
applied around each subsystem. This class of systems is characterized by the 
con5)arison form of tlie subsystem matrices and the lower diagonal form of the 
interconnection matrices. 

Let us consider again the system 
s 


X. = A. 
1 


.X. + y A. *x. + b.u. 
1 1 I ^ ^ 


i = 1, 2, 


(3.Z) 


3 ^ 


■where the n- x n* matrix A- and the n- vector h- are 

XX X X . ^ 
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““0 

1 

... 0 


^ 0 ' 

0 

0 

... 0 


0 






0 

0 

... 1 


« 

0 

^1 

-aj 

2 

... -ai 

il 


1 


and the 




X n . 



matrices 
, p < q 


A.. = fa^^:) 


are such that 


. C3.403 


C3.41) 


■where p = 1, 2, ... , n . and q = 1, 2, ... , n . . 

J 

In order to stabilize system (3.2) characterized by C3.40) and (3*41), 
we apply the local control 

“ -kTx. , C3.5) 

111 ’ 


and get (3.2) as 



C3.42) 


Gain vectors k. are chosen so that each matrix A.- b*kt has a set 1- 
1 111 1 

* 

of distinct real eigenvalues defined by 

'•i “ = ‘““p ! o > 1, > 0, p = 1, 2, . . . , nj. ) 


1 , 2 , ... 


C3.433 


The positive constant a is to be determined, so that the overall system 
C3,2) is stabilized. 

Following the development in Section 3,1, we transform C3.42) into 
s 


X. 

1 


= Vi 


I 




i = 1. 2, 


C3.445 
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■where the transformation C^.IO} is used to get 


A. - b-icT}T. , A.. = 

X X '■ X X X-' X ' xj X X3 3 


C3.45) 


with, in the quT,sidiagonal form 


i i i 

= diag {-atr^, "‘^'^2* *** * ^ “ 


(3.46) 


In this case, the transformation matrix can be factorized as 

■ 


where 


Ui-l 

= diag {1, a, . . . , a } , 


(3.47) 


(3.48) 


and f - is tlie Vandermonde matrix 

X 


T. = 

X 


“CTh 


i 

C-.i3 ^ 


-a. 


_• n.-l 


“i 


C-aJ) " ... C-a^.3 


n.-l 


1 — * 


(3.49) 


For the moment, we consider the free uncoupled subsystems 

x^ " "^i^i ^ X“l, 2, ... >s » 


(3.50) 


Each subsystem (3,50) is stabilized with a degree of exponentisil stability 

-IT. - (3,51) 

where 


IT- - mm a 
1 n P 


(3.52) 
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n. 

1 


Now, we choose again the Liapunov function vj R R.. , 


where 


and 




A^fl. . a. A. = -G. 


Gi = 20^ diag- aa^f ... , acT^_} , L = . 


C3.133 


C3.14) 


C3.S33 


The aggregate system. 


Y < W 


C3,18) 


is formed as in Section 3.1 confuting the elements w- - of the aggregate ma- 
trix W with 


and 


I,. ” xS (aT,a. 0 

M 1] ij-' 




C3.543 


C3.550 


Our ability to stabilize the system depends ultimately on satisfying the 
Sevastyanov-Kotelyanshii conditions (3.21) by the aggregate matrix W = (w- 
defined by 



ij a 






C3.19) 


Since the matrix W has nonnegative off-diagonal elements, it is a well-lcnown 
fact [ 16 i , that the conditions (3.21) are equivalent to the quasidominant 
diagonal property of W , 

> j dilwyl , j = 1, 2, , s C3-563 

ij^j 

where d^’s are positive numbers. Apparently, we can make the matrix W sat- 
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isjfy conditions C^tZO), i£ -we can increase the diagonal elements suffi- 
ciently large idiile keeping the off-diagonal elements bounded. Xhxs is 

exactly the case with the class of systems under consideration. We notice that 
the diagonal elements Ci " 


\j. - 
11 



C3.57] 


depend linearily on the adjustable parameter a . The off-diagonal elements 

Ci ^ j!)» 

, P-585 

are bounded functions of a . To see this, we note that the elements 

o£ matrices RT^A. .R, are either zero for p < q due to C3*41), 
pq 1 ij j 

or they are bounded for p ^ q due to nonpositive powers of a . We have 


lim 

a-H-“ 


R7^A. -R. = D-. , 

1 ij j ij ' 


C3.S93 


where the matrix defined by: d^^ = , when 

d^^ = 0 , when p q . From (3*55} and (3.593, "we define D.. 

pq 

and conclude from 


p = q , and 

= 17^D. .f . 

1 iJ 3 


Ita iyW = C3.6D) 

that the off-diagonal elements w- . are bounded in a . 

Therefore, for the selected class of dynamic systems vs can always choose 
a sufficiently large parameter a , and use local linear feedback control to 
stabilize the systems. From (3.433, see that by increasing the value of 
ct , we move the subsystem eigenvalues away from the origin, thus, increasing 
the degree of exponential stability of each subsystem. This, however, requires 
an increase of the local feedback gains in the course of stabilization. 
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S.4. An IZlus-b^at'iTie Exgrtlple 


Let us illustrate the local stabilization procedure using the following 
exanple: 



C3.61} 


Tlie eigenvalues of the system matrix A corresponding to C3»61), are 

= 1.7244, = 5.1042, X^ = *'1.2633, X4 5 = ~4.2826 ± jl.7755 

C3.62) 


and the system (3.61} is unstable. 

T!ie system (3.61) can be deconposed as 



(3.63a) 


"0 I" 

x„ + 

1 

0 

1 

0 


1 

0 

1 

_-3 -2_ 

z 

L5 

6 


1 

LiJ 


(3.63b) 


The eigenvalues of the subsystem (3.63a) are moved from 
X^ ^ -1.3532, ^2 3 0-1766 ± )1.2028 

to the neiir locations 



applying tlie local control (S.5) and 
^ = (4, 10, 5) . 


(3.64) 


(3.65) 


(3.66) 
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Similarly, the eigenvalues of the siibsystem (3.63b] are changed from 
2 “ ■'1 J 3-4142 , 


(3,67] 


,2 _ 2 _ T ,2 _ 2 _ „ 


(3.68) 


applying the local control (3.5] and 


1^2 - ("1, 1) . 


(3.69] 


Referring to (3.46), we see that in (3.65] and (3.68], the parameter a = 1 , 


We construct the transformation matrices 


R^, R2, t^, T2 for a > 1 as 


10 0 


R, = 0 a 0 , T- = -1 -2 -3 

0 0 14 9 


111 


1 o"i ^ ri 1 

R„ = , T„ “ 

0 a -I -2 


(3.70] 


The numbers h 2 are both set to one. Then, the aggregation matrix of 
(3.18] defined by (3,57] is given as 


(3.71] 


which for a = 1 takes the form 


17.0011 


12.2936 


(3.72] 


^12 " ^*^1 ^?.'^ 2 ^ * ^21 " ^*^2 ■^2 


where 


— -..-..J.,* 


— .1- 


aiid Aj^ 2> specified in C3.63]. 

It is obvious that the matrix If in C3.72) does not satisfy' the in- 
equalities C3.21), 

From (3,63) and (3.59), "we find that 


2 Q“ 

= 0 , °21 
0 0 


4 0 0 

0 ■“ 6 0 


(3.73) 


and for a > 15, , we have 5^2 " 32,55, s 18,98 . Thus, for « - 25 
we have the aggregate matrix 


"-25 32.55" 

W = 

18.98 -25 


(3.75) 


which satisfies the conditions (3.21), and the overall system is stable. The 
corresponding eigenvalues of the overall closed- loop system are 


= -36.0364, X2 3 = -25.9599 ± j3.5219, X^ 3 = -68.5213 ± J6.0474. 

(3.76) 

For the chosen value of ft = 25 , we have the eigenvalue sets and 
defined in (3.43) given as 

h " 


L, = {-aa?, -aah ^ {-25, -50} . 


(3.77) 


The locations of the subsystem eigenvalues specified by ^2 (3*77); 

are achieved by the local state-variable feedback defined by (3,5) and 

= (93748, 6874, 149) 


k2 == (1247, 73) 


(3.78) 


Z8 


The gains in (3.783 are relatively high, which is due to the use of local 
controllers only. The gains can be considerably reduced by applying global 
controllers in the multilevel scheme outlined in Section 3.1 and illustrated 


in. Section 3,2, 
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3,5 AppHaation to 1ST 

In this section, we design a control system for the nonlinear model of 
LSI described in Section 2, by using only the local linear controllers as 
proposed in Section 3,3. This necessitates an application of the results 
obtained by Weissenberger [17] which are concerned with the finite regions 
of stability of large-scale systems rather than their global stability pro- 
perties. 

We notice that the LST model (2.20} belongs to a general class of sys- 
tems described by the equations 

Xi = + a.x]’ A. .X. + b.u. 

fl, 1 = s 




i = 1, 2, ... s; % - 


i+1, i s 


(3.79) 


where A- are constant n, x n. matrices, A. . are n x n. constant ma- 

X IL X XJ] y* ^ 

trices, a- and b. are n. constant vectors. 

1 1 1 

To stabilize the system (3.79), we choose the local control 




' • f 


C3.80} 


so that each uncoupled subsystem 


has a prescribed set of distinct eigenvalues 


(3.81) 


L. = ± ... , a;!; ± , -q 


i 1 


P+1 


, • . . J 


"i’ 


cfg ^ 0; Pj Cl " 2, ... , n^}, i 1, 2, ... , s. (3.82) 


By using the transformation (3,10), the closed-loop system corresponding 


30 


to (3.79) is obtained as 


■*’ * “* ** ~T 

X- = A*x. ^ 

1 XX 


, JL t* ^ 

. + a.Xo )• A..X. , 

I X A X3 j * 

i = 1, 2, 


, s , A == i 


1, i = s 


[i+l, X / s , 


(3.83) 


where A. - T7^(A. - b.kT)T. has the quasidiagonal form (3.12), A.. = 

^ H» IL IL 3* 3L 

tTa. .T. , and a. - TT^a. . 

A xj j ' X XX 

IVe define the interaction function h; T x r“ -> R ^ among the subsystems 
of (3.83) as 




(3.84) 


The interactions h- (x) can be bounded as 


llv^li l^oi .1 r 

j“l 


(3,85) 


on the region 


r = {x 6 R^; x^ < Vq^i^ , i = 1, 2, ... , s} , 


(3.86) 




where are positive yet unspecified constants , and | = (a^a-^) 

?^mCaJ-a..) . 

13 x3"^ 


The aggregate s x s matrix W = (w- •) which corresponds to the system 

ij 

(3.83) and constraints (3.85), is obtained following reference [17], 

W = DW (3. 87) 

where 


D - diag {Vq 2 » Vq 3, ... , ■v'os> 


C3.88) 
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and the s x s :n]atri3c § s ^ defined 


-•1 


13 ij Oii 1 ±j'' 


(3.893 


with 7r^ defined in (3.9). 

Fxoni (3.873) it follows that W satisfies inequalities (3,213 if and 
only if W does. Inequalities (3,21) applied to W determine the constants 
^01' ^02 » * ^Os ^ (3.85). It is possible to calculate these constants 

recursively. To see this, we note that the k-th leading principal k x k 
submatrix can be expressed as 


II 


J ^ 


^ I I 0 
1 

_i ^ 

^k 

1 ""”" 


i 1 

LVk-l ' 


\-i i « 




' i'^-A 

-- 

0 ! 1 


t 


(3.90) 


Therefore, the k-th leading principle minor of W is 


T^rrl 


det Wj^ *- det (Wj^- 


(3.91) 


For the inequalities (3.91) to be satisfied by if , it is necessary and suffi- 
cient that 

~^^kk ®k^k-3^k (3.92) 

From (3.89), we have 

^k “ ^^Ik’ ^2k* » ^sk^ * ®k ” ^^kl» ^k2' *'* * ^ks^ (3.93) 


and from (3,89) and (3.92), we get the constants Vqj^ *as 


.TiTt-I J- ^ “1 


^0* ^ ■’'kCsA-A^ 


S, = 


1 ) k “ s 

k+1, k s 


(3.94) 
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Once the constants are calciilated by C2.94), the region n of 

C3.863 is determined. Now, it remains to imbed a Liapmov function V: -»• 

inside the region n and determine a region of stability fl7] 


n = {x G VC5c) < y} . 


C3.95) 


In (3*95), we choose 

VCic) - ! d-|v.| , 
i=l ^ ^ 


C3.96) 


i^ere d. are positive numbers, and v- = = I 1 • Following [17], 

X aU JL JL 

we calculate the positive constant y it (3*95) using (3.94-] and 


Y = mm i == 1, 2, ... , s 


C3.97] 


T 

where the positive vector d = d 2 » ... , d^} is computed by 

T T — 1 

d = -cV , (3. 98] 

where c is any positive s vector (c > 0) . 

Since ic. = , and | |xJ | < | |t:^| | | [x. ( [ , from [3.96] and [3.97], 

we get finally the region of stability n in the original state space, which 
is 

n = {xGR": f d.||T-^|| ||x.|| < Y> . C3.99J 

i=l ^ ^ 


Now, we consider the nonlinear model of the LST given in Section 2, which 
belongs to the class of systems if ' cribed by (2.79) with 


A. = 

“■o 

11 

1 

, A. . - 

“0 

“1 

■ 3i * "■ 

"0" 

, b. = 

r 0 " 

1 

1 

_0 

oj 

’ 10 

, 

_0 

ij 

1 

tQ. 

1 

’ 1 

j 


i = 1, 2, 3 . 


(3.100] 


-i 


Applying the control law 

~ * i “ 1, 2, 3 


where 


i 2, 3 


C3.1013 


(3. 102] 


and ic^ - , we obtain the closed-loop uncoi 5 )led subsystems C3*81) 


with 


h- hh'^ f 1- , i = 1. 2. 3 ■ 

"•‘^il “ i2 


C3.103] 


The gains 3c^ are chosen so that each subsystem has a set of eigenvalues 


’■o^} > i = 1, 2, 3 . 


C3.104] 


To get the transformed system corresponding to (3,83], we use the trans- 


formation matrix 


1 1 


1 1 
■'^l "'^2 


, i = 1, 2, 3 


(3.105] 


and get 


4 0 1 T- 

-I i, j = 1, 2, 3, i jj £ = 


i+1, i?^3 
1 , i=3 (3,106] 


To con^jute , we choose cr^ = cr^ , Og = ^2 , i = 1, 2, 3 , and cal- 
culate llA^jll = (a^^ + (^ 2 )^ , (a?a^D'^'= l®il^l^l" *^ 2 !^ ^ can mini- 
mize the numbers with respect to the distance p = ^ 2 - between the 

two subsystem eigenvalues. This yields 

lii = |«i| + C®i + P)^] , 


(3.1073 
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j. 


£md i^e get the minimal values | * for L . as 

Xj XJ 


I” “ C4 + > 


[3.108] 


which is obtained for p = , 

Hie corresponding matrix W in [3.87], is 


W = K 


'21 V. 


[3.109} 


From [3.94] and [3,109], we get 


’^Ol ^^02 "^03 


^15 ^21 ^32 


[3.110] 


Choosing Vq^ Vq 2 ‘'^03 ” ^0 ’ using [3.108] and [3.110], we compute 

T 

Vq < 0.584, Selecting Vq = 0.574 , = 10 , and choosing c = [1, 1, 1] , 

we further confute from [3.98] the vector 


d = [4.8, 13.7614, 4.2963)^ . 


[3.111] 


From [3 .97) , we calculate 
Y = 2.4663 , 


[3.112] 


and the region n in the transformed state space as 


n = {^ G K^: 4*8| |x^|] + 13.7614| | + 4.2963||x3!i < 2.4663} . [3.113] 
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In the original space, the stability region fi is finally obtained as 

== {x G R^: 4,8| [x^l I + 13,7614||X2| | + 4.29631 Ix^H < 1.3331} C3.H4} 

where we used | [T^^| | 1.8500 , i = 1, 2, 3 . 

The feedback gains that yield the region J? are coinputed from (3.102^ 

and 

kT = (a^a2» + ^2^^^ = Ci41.4213, 34.1421^'^ , 

i = 1, 2, 3 C3.115) 

as 

k^ = a* 6517, 0.39883^ 

k2 = (10.3303, 2.4939)"^ 

k^ = (10.7056, 2,SM6f . (3.1163 


This completes the design of the LST control system. 


J 


1 


d. 


-I ; . : ; 


I - : J 
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4. OPTIMAL COmOL 


In this section we will describe the application o£ a recently developed 
multilevel optimal control scheme [9] for the decentralized regulation of the 
LST. Such multilevel control schemes are quite efficient in the analysis of 
large-scale systems that may be decomposed into a number of interconnected 
subsystems of smEiller dimensions. Since our model for the LST, described in 
Section 2, is a nonlinear intt^r connected system composed of three linear sub- 
systems describing the motiOiis of the spacecraft along the three axes, gener- 
ation of the necessary control scheme basing the analysis on the subsystems 
is highly desirable in view of the con^lexities involved in the optimization 
of a nonlinear system of a large dimension. In the sequel, we will describe 
the general theory for the multilevel optimal control of interconnected systems, 
whidi will be followed by the specific application to the LST. 

4.1. Frobtem Formulation 

Let us consider a continuous dynamic system described by the differential 
equation 

X " f(x, u) C4.1) 

where x(t) 6 is the state and uCt) 6 is the control fimction of the 
system at time t 6 T . The function f: is continuous on a 

bounded region P ^ R^ and is locally Lipschitzian with respect to x in 0 
so that for every fixed control function uft]) , a unique solution x(t; tg, Xq) 
exists for all initial conditions (t^, Xq) 6 R V and all t 6 T , T being 
an interval [t^, of R . 

We assume that system (4.1) can be decomposed into s interconnected sub- 


= A.x. 
1 1 


B.u. 
1 1 




i = l, 2, 


C4,Z) 


systems 
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n. 

where, ^ is the state of the i-th. subsystem so that 

_ n, n„ n 

R^ = R-^xR^x...xRS; 


m. 

6 R ^ is the control function of the i-th subsystem so that 


m~ ]Tiij m 

R™=R-^xR^x ... xR^ ; 


n. n. n- x n. n 

A. G R ^ ^ and B- 6 R ^ ^ are constant matrices; and h. : R -J- R 

X X 1 

is the ftmction which represents the interconnection of the i-th subsystem in- 
side the overall system. 

The multilevel control scheme [ 9 ] used for the optimization of system 
(4.2) can be developed by considering the control function u^(t) as consist- 
ing of two parts, the local control u^(t) and the global control u|(t) , 


u^(t) = u^Ct) + u|(t) , 


(4.3) 


The local control u^(t) is chosen as a linear control law 

ujct) = C4.43 

to optimize isolated subsystems, and the global control law u?Ct) is chosen 
as a suitable function of the state 


ufm = -KfCx(t)5 (4.55 

to minimize the performance deviation from the optimum due to the presence of 
interconnections among the subsystems . 

With the application of the control (4.3), the equations (4.2) governing 
the system under consideration tahe the form, 


^i 


- A.x. + B.u. 
11 11 


h^Cx) 


B.u? 
1 1 


i = 1, 2, 


(4.6) 


s . 
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Since, as described earlier, the global control functions are assigned 

iXi 

only the task of reducing the effects of interconnections h^fx) the terms 


hg Cx, u|3 = h^CxD + B^u| , i = 1, 2, ... , s , 

may be regarded as the "effective interconnections" among the s 
systems 


C4.7J 


isolated sub 


X. = A.x. 

X 11 


B-uf 
1 1 


i - 1, 2, 


C4.8] 


We shall assume that all s -pairs are completely controllable, 

and that -with each isolated subsystem C4.8) a quadratic perfoimance index 

CO 

JiCtfl, Xj.g, = I { I lx. Ct) 1 1^ + I luj(t) 1 1|_) dt C4.9) 

t ^ ^ 

n. X n. 

is associated. In C4.9) Q. 6 R ^ ^ is a symmetric nonnegative definite 

m. X m. ^ 

matrix and R^6R^ ^ isa symmetric positive defini^te matrix. 

0 

The local control u^(t) in C4.4j can now be chosen to ininimize tlie per- 
formance index ^iO' (4.9)* From linear-quadratic regulator 

0 A 

theory [18], the optimal control (t) is given by 

uf M = -icfxj^W (4.10} 

where 

n. X n. 

In C4.11D> ^ ^ ^ ^ symmetric and is the positive definite solution 

of the algebraic Riccati equation 

P.A. ^ Afp. - PiB.R:^]p, + Qi = 0 . (4.12D 

The optimal cost x^q) = ^j_o» x ^ calculated 





- I. 
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as 


" I l^iol Ip- 


C4.13) 


Furtfiexmorej, under the assun^tson that can he factored as - CX^, 

n- X n. • 

■where C^. G R ^ ^ such that the pair (k^f cp is conqpletely observable, 


each closed-loop subsystem 


.-1t,T. 


X. = (A. ~ B.R. B.P.) X. , i = 1, 2, , s , 

■'1 1 1 1 1 '^ 1 ’ ’ * * ' 


(4.143 


is globally asyn 5 )totically stable. 

oit 

The controls u^, (t) , i = 1, 2, . . . , s , 'will not, in general, be opti- 
mal for the composite system (4X3 and will not result in the optimal cost 


J (tp, Xq3 - I ^0^ 

i“i 


(4.153 


unless the effective interconnection fimctions h (x, u?3 - 0 . "When h (x, 

0* 

u|3 ^ 0 f the controls (t) produce a value of the performance index for 
the composite system given by 


where 


JCtp, Xp3 - J 

1=1 


^iO^ " ^i0» ^i ^ ’ 


(4.163 


(4.163 


It is obvious that 




,n 


(4.173 


.A*, 


and the local control law uv (t3 in (4.103 only be a suboptiraal policy 
for ■the composit!- system (4.63, with an index of suboptimality e > 0 defined 
by the inequality 

*l(tg, Xq 3 ^ (l‘*‘^3 J (^0* ^0^ ^^^0^ ^0^ G T X R , 


(4.183 


The suboptimality index e for the system with, tlis optimal local control, 


^1 K. > ± = 1, 2, 


depends on the size of the effective interactions h u?} 

e^ 1 

measure of the performance deterioration due to these. 

We can now give a formal definition of this concept. 


. , s C4.19} 
and hence is a 


Definition, The system (4.19) with the optimat local control law (4.10) is 
said to he suhoptimal with the index e if there exists a number e > 0 for 
which inequality (4.18) is satisfied. 


As described earlier, the suboptimality index e is a function of the 

interactions h (x, uf) and the following problem is of interest r 
®i ^ 

"Brohlem 1, Establish conditions on \i (x, u?) to guarantee a prescribed 

H ^ 

value of the suhoptimality index £ . 


It is important to note that in Problem 1, the rate of the global control 
function u|(t) is ignored as it is taken together with the existing inter- 
connections h^(x) in the system to yield the effective interconnection func- 
tion h (x, uf) . However, as we shall see later, the solution to Problem 1 
indicates a method of choosing the global control u|(t) so as to reduce the 

size of h (x, u?) and, hence, minimize the suhoptimality index e . In 
®i ^ 

other words, we consider the index e = e[| |h^(x, u®) 1 1] where h^: x 

+ r" is hg = [h^ , hi , ... , ]’’ iv.d is uB = 

p T T ^ ^ ® 

... , (u°) ] and solve the following; 


,m 


Frohlem 2. Find a control law of the form (4.5) or equivalently-^ 


u°rtj = 


for which. 


C4.2Q} 


1 ! 1 ] J 1 


e = inf e{[||h (x, u^)||]}VxG'P 

w 

u^C^) 

is attained. 

We will now provide the solutions to the above problems. 


.J i 
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(4.21) 


4,2, Multilevel Optimisation 

A solution to Problem 1 may be obtained by using the classical Hamilton- 
Jacobi theory. Since in our optimization procedure, we chose the local con- 
trol laws (4,10) to optimize the decoupled subsystems, the optimal indices 
satisfy the corresponding tiamilton- Jacobi equations. When the subsystems are 
interconnected, the equations are not satisfied by the respective performance 
indices and the overall system is not optimal. However, a raajorization pro- 
cedure is possible to provide an estimate of the performance deviation from 
the optimum due to the interactions. 

Now, we provide a solution to Problem 1 by the folloidng; 


Theorem. 4.1. Let there exist nonnegative numbers . such that the function 
h (x, uf) in (4,19) satisfy the constraints 

6 « X 

1 

ll\ “&II - ^ SinIKII . VxER’' , Vi = 1, 2, ... , s C4.22) 
i j“l *' 

and 





VK) 


C4.23) 


s s 


where K = I I ?. ■ ,P= diag{P_ , P„, 
i=l j-1 ^ ^ 


, P^} , W- diag{W^, W 2 , ... , W^} , 


P^ being defined by ('4,12) and ~ 

are the maximm and minimum eigenvalues of P and W respectively. Then the 
composite system (4.19) is 

(i) suboptimal with index e 
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and 


(ii) globally asymptotically stable. 


Proof, Since the decoupled subsystems (4.14) are optimal, the functions 

V. (x.) = I [x. 1 Ip , i = 1, 2, ... , s , satisfy individually the Hamilton- 
11 1 

Jacobi equations 


[grad v^(x^)]'^[CA^- 




n. 

1 


Vx^ G R , i = 1, 2, . . . , 


(4.24) 


Now, the time -derivative v. (x.) can be calculated along the trajector- 


1' 1 


ies composite system (4,19) as 


ViCJp = [grad V C%)]" {(A.- B.R'^E^Pp 5^ + (S, u|)> 

i 

, „ -..T »T ~TiT 

where x = [x^, X 2 , , x^] . 

Substitution of (4.25) in (4.24) and rearrangement of terms gives 

11x^1 1 [grad v^(x-)]'^ hg (x, u|) 

i i 

-e||x.|1^, Vx S , 1 = 1, 2, ... , s , (4.26) 

i 

where the simplification \ \x^\ |? + | |K^ x^| = | [x^| with " Qi 

P . B . R. . is made . 

Ill 11 

s 

Denoting v(x) = V v- (x.) and summing the s- equations in (4.26) we getj 

• >1 IX 
1=1 

I lx| Iw = " d+e) v(x) + (l+e)[grad v(x)]'^ h (x, u®) 


(4.25) 


- e||x||2, VxGR’^ . 


(4.27) 
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NoWj integrating C4.27) from t^ to » we obtain 


+ Cl+e) 


{[grad vCx)]"^ h„Cx, u^) - tI^ I l==l lw> 


l-»-E 


I 


^0 


where J and J are defined in (4.15) and (4.16]). 

It is now simple to observe from (4.18) and (4.28) that the system is 
suboptimal with index e if 


{[grad v(x)]'^ hg(x, u^) - ^ | jx] |^} dt > 0 , 


Vx e . 


(4.29) 


For further simplification of (4.29) we note that 


v(x) = I " J ll^llp ^ ll^llp . 

i=l ^ ^ i=l ^ ^i ^ 


(4.30) 


s 

Also, since ||h (x, u|)|| < J g. . | [x. | |,Vx 6 R we have the inequality 

^ ^ j-1 J 


I u®)| I <_ cl |xl [ , Vx e R^ 


(4.31) 


s s 

where ? = I I ■ 
i=l j=l 


Using (4.30) and (4.31) it can be easily shown that a sufficient condi- 
tion for the inequality (4.29) to hold is 

ZSKl 1*1 1 i xl? 1 1^1 Iw * ® 

which, however, is in^jlied by the main inequality (4.23) of the Theorem. 

To con^lete the proof of the Theorem, we demonstrate the global asympto- 
tic stability of the system (4.19) by using the function v(x) = 1 [xHp 


as a 
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] I._ 1 =1 


! 


Liapunov function. Note that v(x3 is positive definite since P is a diag- 
onal matrix formed from the positive definite solutions of the s Riccati 
equations (4.123. Further, the time-derivative of v(x3 along the solutions 
of (4.193 gives 

V(x3 = - I |x| + 2x'^Phg(x, u®3 <0 Vx G , (4.33} 

from (4.31} and (4.233. This con^letes the proof of the Theorem. 

It is important to note that the above theorem provides an eiqplicit al- 
gebraic constraint on the interactions that is easy to check. Inequality 
(4.23) involves calculation of eigenvalues of block-diagonal matrices P and 
W , and since = m^ \j^(W) = min j the calcula- 

tion can be carried out on the subsystem level. 

In the context of the above Theorem, it is of interest to consider Prob- 
lem 2 of determining the global control u®(t) so as to minimize the subop- 
timality index e . From (4.233 and (4.31), it is evident that e is a non- 
decreasing function of j [h (x, u^) | [ and hence, Problem 2 reduces to one of 
choosing u^(t3 to minimize ||h (x, u^3ll . This function minimization prob- 

w 

lem is particularly simple in the present case since, from (4.73 

h_(x, u®3 = h(x3 + Bu® (4.343 

which, on using the control law (4.53 reduces to 

h^(x, u^3 = h(x3 - BK®(x3 . (4.353 

A perfect neutralization of the effects of interconnections occurs if a 
choice of K®(x3 results in 

BK®(x) = -h(x3 C4.363 

and, in this case, e 0 . In the special case, when E is square and non- 


45 

singular, the explicit expression for is available as 

kS(x 3 = -'B‘^h(x3 , C4.37) 

In general, a perfect neutralization of the interaction effects mentioned 
above, is not possible and one may attempt to minimize | |hCx) - EK^(x)| | by 
the proper choice of K^(x} in order to solve Problem 2. This is adinittedly 
a conplex minimization problem and a general solution is difficult to obtain. 
However, in the particular case of linear interconnections, the problem can 
be simplified and an elegant solution can be provided. This is, we assume 

h(t, x) = Hx (4,38) 

where H 6 , In this case, the global control can also be chosen as a 

linear law 

K^(x) = K^x (4.39) 


where 6 . With (52) and (53) , Problem 2 sinqjlifies to: 

Prohlem s\ Choose the matrix K® such that inf 1 1 (H-BK^)x| | is achieved 

n KS 

for all X S r" . 

Remembering that ( [ (H-BK®)x| j £ [ jH-BK^l ( j [x( | holds for all x 6 R^ , 
Problem 2’ actually reduces to finding min j |H-BK®] | . When rank B = m , the 
solution to this latter problem is well-known and is given by 


I K® == (B^)‘^ B% 

i 


(4.40) 


T “1 T 

where (B B) B is the Moore-Penrose generalized inverse of B [15] . It 
is interesting to note that in the particular case when 


Rank [B 1 H] = Rank B 


(4,41) 
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the choice (4.40) leads to a perfect neutralization, of interaction effects 
and e = 0 . 


4.5. An ItluGtrative Example 

For the purpose of illustrating the multilevel control scheme presented 
here, let us consider the following exanple. 

The system is described by 

X = Ax + Bu • C^.42) 

where 



-5 

6 

0 

-0.0095“ 



"l 

o" 

A = 

4 

-4 

0.003 

0 

and 

B = 

0 

0 


-0.00332 

0 

-3 

1 



0 

0 


- 0 

0.00995 

8 

-2 



_0 

1_ 


and is required to be optimized with respect to the performance index 


J = 


{||xir •»* iiuir} dt , 


(4.43) 


In this particular case, the dimension of the system 0^ = 4) is small 
and hence the problem is amenable for a direct analysis and the required con- 
trol can be obtained from solving the associated Riccati equation (of fourth 
order). However, as our interest here is to provide an illustration of the 

A 

decentralized optimal control scheme , let us consider system (4-42) as being 


* 

Besides the advantage of permitting an analysis based on the subsystems of 
small orders, the decentralized control scheme presented results in important 
connectivity properties of the system* The suboptimality and stability of 
the system remain invariant under structural perturbations caused by the on- 
off participation of the parts df the system. This property, however, does 
not result when direct optimization of the system is carried out tlb] . 
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composed o£ two subsystems 

C4.44} 


C4.4S) 


with the interconnection matrix 

“ 0 0 

0 0 

^ "" -0.00332 0 

0 0.00995 


0 

0.003 

0 

0 


-0.0095 

0 

0 

0 


C4.46} 


Xi = 


and 




15 

6~ 

X- + 

1 

[ 

1 — 

-4_ 

1 

Loj 

-3 

1~ 

x„ + 

1 

^8 

-2_ 

£1 

— r 

t— i 


Ui 


u. 


B)’’ splitting the control functions and into a local conponent 
and a global conponent, the decoupled subsystems C4.44D and (4.45) can be op- 
tiinized with respect to the performance indices 


{||x^||^ + llu^ll^} dt and 


t 


0 



(4.47) 


The solutions of the associated Riccati equations can be obtained as 


“1,1910 

1.5411 

and 

P„ = 

5.5591 

2.3746“ 

_1.S411 

2,1397_ 



_2.3746 

1.1224_ 


and the local control laws are, 

u^ = -[1.1910 1,5411] x^ 

U2 = -[2.3746 1.1224] X 2 (4.48) 

In the absence of the global controls, the functions (4.48) will only 
be suboptimal policies for the overall system (4.42) with the index of sub- 


optimalily e given by 


liHlI 



2 max{Aj^(P^} f ^^^2^ ^ 


C4.49:i 


where W. = Q. + P.B.r7^bTp. , i - 1, 2 , are 

1 1 1 X 1 0. ’ ’ * 


2,4185 

1.8354'" 

and 


-6.6386 

2.6651 

_1.8354 

3.3748_ 


2 

_2.6651 

2.2597_ 


Inequality (4.49) is satisfied with e == 2 and hence the performance 
degradation from the optimum is 200% , 

In order to improve the performance, we now use the global controls u| 
and u| given by 

b'^Hx; (4.50} 

where u® = 


U‘ 


g 


u; 


g 


. (4,50} can be sii^lified to yield 


uf=- 

E 0 

-0,0095 } 


uf = - 

[ 0 

0.00995] x^ . 

(4.51} 


Tne effective interaction matrix fl with the application of the global con- 
trol is 

H - [ I - B(b'^B}"^ b"^] H 

^ 0 0 0 0 

_ 0 0 0.003 0 

TOO 
0 0 0 


-0.00332 

0 


(4.52) 
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and the suboptimality inequality (4.49) can now be satisfied with a value 
e == 0.2 . Hence the performance degradation is reduced from the original 
200% to only 20%, thus illustrating the effectiveness of the global controls. 


4.4, Application to LST 


The results developed in the earlier parts of this section may be di- 
rectly used for the multilevel optimization of the LST, As described in Sec- 
tion 2, the model for the LST is a set of three interconnected subsystems, 
described by (2.20) as, 

^i ” ^i^i ^ ^i^i ^ » i = 1, 2, 3 , (4.53) 



b 1 


"o" 


0 n 


2 

where x. G k , A.= 

0 o_ 

. b. = 


and hj^ (x) - 

i 

] 

7 

X X 


J“2^32^2 J 



h^Cx) = 


■“ 2 ^ 32^12 


, hjCx) = 


■“ 3 ’^ 2^22 


, X being the composite state- 


T T T T T 

V0ctoi? 2C “ r ^ ~ ^ ^ 


Following our multilevel control policy, we split each of the control 
functions u. into a local component u^ and a global con^onent u| and 
optimize the decoupled subsystems 


X. = A.x. + b.u. 


1 1 


1 1 


i = 1, 2, 3 , 


(4.54) 


with respect to the performance indices 



i = 1, 2, 3 , 


(4.55) 


obtained with the choice ^i ” ^ Vi = 1, 2, 3 . The solu- 

tion of this linear-quadratic optimal control problem is simple and involves 


the solution of the associated Riccati equations, 


aTp. + P.A. - P*h*bTp* + Q. = 0 , i = 1, 2, 3 
11 11 1111 u 


(4. 56] 


With the specified structure of A^ and , the solution of (4.56] can be 
obtained as, 


P. = 
1 




Cl* 






2 .H 




C4.57) 


and the local optunal controls are 


"i “ ■'''iVi “’Cl Cl + 57)^] , i = 1, 2, 3 . 


(4.58] 


However, in the absence of a suitable choice of the global control func- 
tions u| , (4.58] will only be suboptimal for the composite system (4.53), 
with the index of suboptmality e determined by the size of the effective 
interconnections , 


h (x, u?) = h- (x] 

w* X X 

1 


+ b.u? 


1 1 


i = 1, 2, 3 . 


(4.59] 


(4.59) can be simplified to yield 


^ej^Cx, u|) 

\Cx, <|) 

\Cx. ufD 


0 


2“l^2Z^32 

0 

"“ 2 ^ 2^32 

0 

'“3^2i^22 


"^1^1 


*®2i 




C4.6tt3 


It is now simple to observe that the choice o£ 


“fW = >^2*32 

^ x^2='32 


C4.61) 


will mahe the functions h (x, uf) = 0 and hence e = 0 , thus resulting in 

0 • X 

1 

no degradation of the performance from the optimum. 

It is of interest to evaluate the control functions for a representation 
set of values of the parameters of the LST. For the values of the inertia 
components* = 14656 Kg^ , 1^ = 91772 Kg^ and I^ = 95027 Kg^ and typical 

= = 12.57 X 10^ N-m/rad. , the values of 

a. ,6. , i = 1, 2, 3 can be calculated as 

= 0.2221 = 85.62 

= -0.08754 and B 2 “ 13.69 

a- = 0.8112 = 13.21 

o o 

Hence, the control components can be evaluated from (4.58) and C4.61) as. 


reaction wheel constants 




= K. 


uj = -[1 

1 , 012 ] 

^11 

_^ 2 _ 

“2 = -[1 

1 . 061 ] 

^21 

_^ 22 _ 

“3 = 

1.07 ] 

“^sr 

_^32 1 


C4.62) 


These values correspond to the on-orbit configuration of the LST with extended 
light shield and solar wings, with the corresponding mass of the body totalling 
9380 Kg [11] . 
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and 


« 0.0026 ^2Z^32 
u| = -0.064 3 Cj^2^32 
uf ” 0.0613 3C^2^22 


C4.633 


TMs coinpletes the multilevel optimization of the LST control systsm. 
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.1 


I 


5. CONCLUSIONS 

A multilevel scheme was proposed for control of Large Space Telescope 
modeled by a three-axis-six-order nonlinear equation. Local controllers were 
used on the subsystem level to stabilize motions corresponding to the three 
axes. Global controllers were applied to reduce (and sometimes nullify) the 
interactions among the subsystems. A multilevel optimization method was de- 
veloped whereby local quadratic optimizations were performed on the subsystem 
level, and global control was again used to reduce Cnullify) the effect of 
interactions . 

The proposed multilevel stabilization and optimization methods are pre- 
sented as general too^-s for design and then used in the design o^ the LST 
Control System. Furthermore, the methods are entirely computerized (Appen- 
dices A.l and 2), so that they can accommodate higher order LST models with 
both conceptual and numerical advantages over the standard straightforward 
design techniques. 
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APPENDIX 

COMPUTER APPLICATION 


A.l. Stabilization Program 
A. 2 . Optijnization Program 


A.l 
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A.l« Stabilization Program 

The entire stabilization method is computerized. In this section we 
present the conputer program for tlie stabilization of a class of large-scale 
systems by local state feedback, according to Section 3.3. The program is 
written in FORTRAN IV, for 1IP2100 conputer. It is, basically, an interactive 
user oriented program. 

Designers can enter the program from VDD (Visual Display Unit) and thus 
freely alter the course of confutation, according to the nature of the prob- 
lem. Program accepts input data from a logical unit that has to be previous- 
ly assigned. As a result of computation, it prints out stabilizing parameter 
a , corresponding aggregation matrix, stabilizing set of subsystem eigenvalues, 
and enables the designer to reenter the program with so confuted new set. The 
program finally prints out the corresponding subsystem feedback gains. The 
name of the main program is PPl. Its function is to coordinate the sequence - 
of actions during the course of the stabilization and to enable the designer 
to access the program at various points during its operation. The program 
PPl calls subprograms, DECP, PPL, TRF, AGR and MINV. The processing of vari- 
ables beticeen the main program and subroutines is realized via CCWMON block. 

Program PPl 
Purpose! 

Local stabilization of a class of large -'Scale linear systems. 

Description of input parameters: 

A - N by N system matrix, 

B - N by M input matrix, 

II - one dimensional integer array. It stores dimensions of each 
subsystem. The other parameters are working variables. 


A. 2 


User has to specify integers N and M and nmber of subsystems IS, 
into vMch system matrix A , and input matrix B are decomposed. During 
the course of stabilization, use'" has to enter the program ■with subsystem 
eigenvalues, and specify an increment delta by idiich a is increased dur- 
ing the process of iterations. 

At the very beginning of the program, the user has to assign input-out- 
put units. Also, during the operation of the program, user communicates with 
■the program by speci:fying commands, by ■^^^lich the sequence of calculations is 
controlled. These commands are in the "question- answer” form. For exangsle, 
program prints out the question: 

"DO YOU WANT TO CONTINUE, YES OR NO". The user then types either 'YES" 
or "NO" accordingly. Other commands are self explanatory, and are not going 
to be discussed here. 

Subroutine DEGP 
Purpose: 

Decomposes system matrix A and input matrix B into subsystems . The 
product of the decomposition is stored in A2 and B2. 

Usage: 

CALL DECP CIS, M, N) 

Description of parameters i 

IS - Number of subsystems. 

M - Number of inputs. 

N - Order of the overall system. 


The following parameters are passed via CdM)N block as; 
COMMON A, B, A2, B2, II 
A - N by N system matrix. 
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B - N by M input matrix. 

A2 - Three dimensional array •which, contains the product of 
the decomposition of the matrix A . 

B2 - Three dimensional array "^diich contains the product of 
the decomposition of the matrix B . 

II - One dimerisional integer array vMch contains the dimen- 
sions of each subsystem. 

Subroutines required: None, 

Subroutine PPL 
Ptoppose: 

Pole shifting using state feedback. 

Usage: 

CALL PPL CN, m 
description of parameters: 

N - Order of -the system. 

BV - Integer for the output logical unit. 

The follo^ving parameters are passed iria C0^MDN block as: 
CQ^^vI0N Al, Bl, B2, II A, Ql, Q, 111, B, Rl, R2, D, SK 
A “ N by N system matrix 
B " N-'th' dimensional input •vector 
Rl “ One dimensional array which contains real parts of 
eigenvalues of matrix A . 

R2 - One dimensional array wiiich contains imaginary parts 
of eigenvalues of matrix A . 

D - N-th dimensional gain vector. 


All o-ther parameters are trarking variables, which are placed in COM^fON 
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block in order to make it consistent witB. the CQNWON block of the main pro- 
gram PPl. 

SubrouHnea required: ALAM, DISP, KBAR, 

As a result of tlie pole shifting, the subroutine passes back matrix A 
of the closed loop system, the gain vector D , and the nevf eigenvalues. The 
subroutine itself is inritten as* a user-oriented interactive program^ The user 
enters the desired eigenvalues from VDU, The Coplands for controlling a se- 
quence of computations, are self explanatory. 

Subroutine TRP 
Purpose: 

Transforms subsystems by similarity transformation. 

Usage: 

CALL TRF (IS) 

Uesaription of parameters: 

IS - Number of subsystems. 

The following parameters are passed via COMMON block as: 

C0^M)N A, B, A2, BZ, II, Al, Ql, Q, 111, B3 
A2 - Three dimensional array. It contains the product of 
the decomposition of the system matrix' A , 

B2 - Three dimensional array. It contains the product of 
the deccmposition of the input matrix B . 

II - One dimensional integer array that contains dimensions 
of each subsystem, 

Q “ Three dimensional array that contains transfoimation 
matrices. 


All other parameters are working variables. 
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Subroutines required; MlfsV 

The product o£ transformation is in A2 and B2. The array Q is un- 
changed. 

Subroutine AGR 
Purpose; 

Forms an aggregate matrix. 

Usage; 

CALL AGR CIS) 

Description of parameters: 

IS - Number of subsystems. 

The following para:neters are passed via CCtM)N blocJc as: 

COMMON A, B, A2, B2, II, A4, A3, Q, 111, B3, Rl, R2, KB, SK 
A2 - Three dimensional array idiich contains the product of 
the decomposition of the matrix A . 

B2 - Three dimensional array which contains the product of 
the decomposition of the matrix B , 

II - Integer array that contains the dimensions of each sub- 
system. 

Rl - One dimensional array that contains real parts of sub- 
system eigenvalues, 

A3 - Matrix that contains the aggregate model. 

All other parameters are working variables. 

Subroutine required: ALAM, BIGl, SMALl, 

Subroutine KBAR 
Purpose; 

Computes gain vector for state feedback control. 


Usage: 

CALL KRAR CA, N, Z, IZ, D, B] 

Description of parameters: 

A - N By N system, matrix. 

N *• Dimension of the system. 

Z “ One dimensional array that contains the desired cliarac-^ 
teristic polynomial, 

IZ - Its dimension. 

D - One dimensional array that contains resultant gain vector. 
B ” Input vector. 

Subroutine required: COEFl, SCALU, VECPR, MIW. 

Method: 

Described in reference [13] . 

Subroutine ALAM 
Purpose: 

Calculates eigenvalues of general N by N matrix. 

Usage: 

CALL ALAM (A, N, D» COF, Rl, R2) 

Description of parameters: 

A - N by N system matrix. 

N - Dimension of the system. 

D " N+1 dimensional working vector. 

COF - N+1 dimensional working vector. 

Rl " One dimensional array of real parts of eigenvalues 
of matrix A . 

R2 - One dimensional array of imaginary parts of eigen- 
values of matrix A . 

Subroutines required: COEFl, POLRT 
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Method: 

Confutes coefficients of characteristic polynomial, and calculates its 
zeros. 

Subroutine CQBFl N, D) 

Piirposs.* 

Calculates coefficients of the characteristic polynomial of 
matrix A . 

V&age: 

CALL COEFl CA, N, D) 

BesGription of pcwmeters: 

A - N by N system matrix, 

N " Dimensions of the system, 

D - One dimensional array of coefficients of characteristic 
polynomial. 

Subroutines required: UNITl, PRODl, TRACI, SCMl, ADDl 
Method: 

Uses Souriau-Frame-Faddeev algorithm. 

Subroutine DISP 
Purpose: 

Form polynomial from its zeros. 

Usage: 

CALL DISP CRl, R2, N, Z) 

Descriptions of parameters: 

R1 - One dimensional array of real parts of roots of a given 
polyncmial. 

R2 - One dimensional array of imaginary parts of roots of a 
given polynomial. 
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N - Order o£ a polynomial. 

Z - One dimensional array* that contains computed coefficients 
of the polynomial. 

Subroutine required: FMPY 

All other subroutines, listed in the Appendix are self explanatory and 
are not going to be e3q)lained here. Subroutines PMFY, POIRT and MINV are 
IBM-SSP subroutines. 


:pRiFrTN4.99 

FTM'fiL 

PROGRAM PPX 
C 

C INTERACTIVE PROGRAM TO STABILIZE A CLASS OF LARGE SCALE SYSTEMS 
C 

c . 

c 

c 

C AUTHOR i MIROSLAV B. VUKCEVlCi DATE /7/20/197S/ 

C 

c 

REAL KB 
INTEGER AlO 

integer dd 

DIMENSION AtlOtlO}«B(lOt31tA2(9«IOilO}t 
1 B2O,10«3] 1 1I tSlfAKtSiS] tQI (5t5) ,Q(3t5t5} ,111 (5*5) ,B3(5) ,R1 (6) , 

1 R2(G] fKBCE) , SKIS, 101 ,A1 t3,10},Tl C6,6) •Y2{&,6) ,MHt25) «LLt251 ,CtE5) 
CQKKON A,B.Aa,32tII,AK,Ql,0iUl,B3,Rl,R2iKB,SK,Al 
DATA D0/2HYE/ 
kR1TE(1,100U) 

1000 FORMAT {10A»"STABIL1ZATI0N OF A CLASS OF LARGE SCALE SYSTEMS”) 
WKITEU*1AS) 

IA5 FORMAT! lXf"ASSIGN LOGICAL UNITS"/1X,«212") 

C 

C READ THE DATA 
C 

READlI«lllllRD«lN 

HHITEttH,U6) 

1A6 FORMAT tlX»”SPECIFY ORDER OF THE SYSTEM AND NUMBER OF INPUTS” 
tlX,''2I2”J 
REAOtlRUaillNiH 
111 format ISI2) 

KRITEtIij,150) 

ISO format t IX, ”EHTER SYSTEM MATRIX”) 

DO a K=I,N 
WRITEtI*.I6I) 

IS! FDRMATClXt”, ,”,n, 

I ,”) 

8 REAUtIRD*100)CACK*JI*j£l«N) 
iiRITEtlNflSZ) 

152 FORMAT (IX, "ENTER INPUT MATRIX”) 

00 9 Ksl.N 
WRITE CtW, 151) 

9 BEAOrlRD.100)lB(lt,J),J=l,H) 

* 100 FORMAT (SFIO. 01 

C 

C WRITE THE DATA 
C 

WHITE tlK,530))r,H 

530 format tlX»”ORDER OF SYSTEM * "»I2, "NUMBER OF INPUTS ■ ”il?) 
hRUE(Ii,«53E] 

532 format (ISX, "SYSTEM MATRIX”) 

GO 3T1 K=1,N 

371 W«IT£CJW,105) (AtK,J)»J»l,N) 

WRIT£(IW.B37) 

SJ7 FORMAT (ISX, "INPUT MATRIX”) 



PAGE IS 


oAr.f noo) 



UUU uuu ooo OOO OOO tjoti 


«. DO 372 K«liN 

372 tfRITEllH.105) CBtlC»3)»J=l»Hl 
105 rORHAT(lXfSFl<*.6I 
- WRITE CIW, 147) 

147 FORMAT tlX*«SPECIFT NUMBER OF 5UBSYSTEM5'«/1X.”12'') 
READfIRDflllUS 
W«ITEtIVi»544nS 

544 FORMAT nXt ''NUMBER OF SUBSYSTEMS - "tIZJ 

WRITEClWtl49) 

149 FORHAT(lXt"SPECIFY ORDER OF EACH SUBSYSTEH"/lXf"5IE») 
READURDf nil {Il(K)tK=liI5) 

DO BOO K^nis 
800 WRITE (IW»554SH<tII{K) 

5545 FORMAT I lX»t'0HDER OF SUBSYSTEM". I2i"="tI2) 

START DECOMPOSITION 

215 continue 

DECOMPOSE SYSTEM INTO SUBSYSTEMS 

CALL DECPtlS.H.Nl 

START STABILIZATION 

DO 10 Ksl.IS 
Ll=IItK) 
lP=tK.*lJ 4lS.lt 
00 20 L*1.L1 
00 20 Jn.Ll 
20 AKtL*J)=A2tlPrL.J) 

DO 30 L-ltLl 
30 03(L)=fa2(KiL.Kl 

LOCATE POLES OF EACH SUBSYSTEM 

CALL l>PL(LliIW) 

DO 15 L=l.Ll 
15 TKLtKIsRltL] 

10 CONTINUE 

wRITEdWtSOOl 

500 format t IX. "DO YOU WANT TO CONTINUE. YES OR NO") 
READd.SOllAlO 

501 FORMAT tA2J 
IFIAID.NE.DOIGO TO 2G6 

START ITERATION FOR ALFA PARAMETER 

wRlTEdW.600} 

600 FORMAT dX.«SPECIFY iNCSEHs nELT"/lX."F10.0") 
REAOd.EDOlOELT 
200 FCRMATIFIO.OJ 
ALFal, 

211 CONTINUE 

FORM VANDERMONDE MATRIX 


00 112 K=l.IS 
Ll=IItNI 


PA'JE 1)002 



1 

! 
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noo non ooo 


i>ARF nuoj 


00 112 J>1>U 
SK«JfK)-Tl tJ*K) 

112 SKtJ.Kl*ALF«SK(J.K) 

DO IZ K*l*15 
Ll*IItKl 
□0 OB L>lfLl 
00 B8 JsItLI 

BB 0C!«*L«JJ=5K(JiKl**tL-l) 

12 C0><TINUE 

FORM AGGREGATE MODEL 

call DECPUSfHtM) 

CALL THF{1S) 
call AGHtlS) 

HRITEtlWiB31] 

531 FORMAT aXiXAGGREGATlON MATRIX") 

DO 281 K=U1S 

Zai WRITE t6t 105) tQltiCtJ>tJ=lt IS) 

ALF=ALF*DELT 

CHECK THE SIGN OF THE K-TK MINOR 
1Z*1 

DO 2B2 K^ltlS. 

IF=K 

DO 283 I=liIF ^ - 

DO 283 J=»fIF 
IP = tl-l)S'K*J 
283 C(lPj=QHJ*l) 

CALL MlNVtC»lF.D,LL»HK) 

IZ*-I2 
■ D=IZ*D 

wRITEt6t270)D 
270 format tlXfFlA,6) 
lF<0)211»EllfEB2 

PRINT CORRESPONDING EIGENVALUES 

2B2 CONTINUE 

DO 11 K=lilS 

LlslltK) . . - - 

WRITEClWrTQOJK 

700 FORMAT HX*"EIGENVALUES OF SUBsYSTEM'MZ) 

DO 18 JsliLl “ - . 

18 ■RlTEC6»701M»SKtJ.K) . / 

701 FORMATtlXt"LAMaOA{«»lZ»"J ■ "iFl5.6) 

11 CONTINUE - * - - 

GO TO 215 
268 CONTINUE 
END 
ENOS 
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*■ . 

0161 


SUBROUTINE BIGKRltLiN) 


0162 

C 

. 


0163 

c 

PROGRAM TO CALCULATE THE 1 


0164 

c 



0165 


DIMENSION Rlt6) 


0166 


K*1 

•• ■ 

0167 



- 

0163 


.10 1F(HUK>.GE.R1(J1}G0 TO 1 


0169 


K“J 

• 

■ 0170 


8 J*J*1 

-■ 

oin 


IftJ.LE.NlGO TO 10 

*■ 

0172 


L*lt 


0173 


RETURN 


0174 


end 


•* NO. 1 

E3R0H5*r PROGRAM * 00055 




> 

M 

to 


% 
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0175 SU5R0UT1RE TRANl (Ql i AG tLI iLS) 

0176 C 

0177 C PROGRAM TO TfiANSPOZE A MATRIX 

0178 C ■ 

0179 DIHENSIDK OKSii tAGtS.S) 

0180 DO B KsltUl 

OlBl DO B J«1,LS 

0182 8 A8(J*K)«Q1 (Kf J) 

0183 RETURN 

0184 END 


•* NO ERRORS** PROGRAM * 00056 COMMON » 00000 
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0105 


SUBROUTINE SMALl (RltLiN) 


OiSb 

C 



0187 

c 

PROSRAM TO CALCULATE THE SMALLEST 

element of an array 

0188 

c 



0189 


dimension Rl{b) 


0190 


K*1 


0191 


jan*l 


019Z 


10 lFtRl(K),LEiRlU))GO TO 8 


Ol93 


9 K*J 



0194 


a J»J*1 


0195 


IFtJ.LE.NlGO TO 10 


019b 


L*K 


0197 


RETURN 


0198 


END 



•• SO ERRORS** 


PROGRAM > 00055 


COKHON a 00000 


019 V 
OZOO 
D201 
0ZP2 

0203 

0204 

0205 
0204 

0207 

0208 

0209 

0210 
0211 
0212 

0213 

0214 

0215 

0216 

0217 

0218 

0219 

0220 
OZJl 
0222 

0223 

0224 

0225 

0226 

0227 

0228 

0229 

0230 

0231 

0232 

0233 

0234 

0235 

0236 

0237 

0238 
02J9 

0240 

0241 

0242 

0243 

0244 

0245 
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SUBROUTINE K8AR(A»NfZ* IZtDiB) 

C 

C PROCHAH TO CALCULATE GAIN VECTOR 
C 

OINENStON A(5t5) fZU) •D<6) iB 15) tQlSfS) «R(5) »P1 (23) tLLtS) tHN(5) 
lfCl5) 

CALL COEFKAtNfO) 

DO 8 KeIiN 

8 Z(K>■Dt^)-Z(K) 

□0 9 I>1.N 

9 a(l»N)EB(l) 

L*N-l 

DO 99 Si»liL 
L1=N-K«1 
L2=N-K 
F»D(LU 
DO 7 IeI.N 
7 aiI)®Qn*Nl 

CALL SCALU(BiNiF) 

DO 77 I*1,N 
77 ftUf«Q(l.Lll 

CALL VECPB(AtRtNiC) 

DO I’O Ixl,N / 

10 QtlfL2]=C(I)*8(l) 

99 CONTINUE 

DO 200 «*ltN 
DO 200 JsliN 
Ia(K-l)*N*J 

200 Pl(IJ=aUtK) _ 

DO lOS KaltN 

lOS B(K)=U(KtN) 

CALL HINVtPl.N»SD.LL»MM) 

IFtSU)S»4»5 

4 WRITE (Ski 02) 

102 FORMAT aXtl3HSINGULAR CASE) 

5 CONTINUE ■ " ” . 

DO 201 K>liN 

DO 201 J«ltN _ 

I»tK-l)*N*J 

201 Q(J*K)eP 1(1J 
DO 91 2»ltN 
S»fl« 

DO 92 Jsl«N 
92 S>SoZ{J)»Q{JiK) 

91 DtX)>S 
RETURN 
END 



*• NO ERRORS *4 
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COMMON X DOOOO 
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0246 

0247 

0248 

0249 
0260 
0261 

0252 

0253 

0254 

0255 

0256 

0257 
0253 

0259 

0260 
0261 
0262 

0263 

0264 

0265 
0c66 

0267 

0263 

0269 

0270 

0271 

0272 

0273 ' 

0274 
02 f 5 

0276 

0277 
0273 
0279 


SUBROUTINE DECPtIS«HiN) 

C 

C PROGRAM TO DECOMPOSE SYSTEM INTO SUBSYSTEMS 
C 

DIMENSION A(10tlO)iBUOO}iA2{9»iatlO]fB2(3flOO) iIKB) 

COMMON AtB«A2it)2*ll 

IL»1 

L1«0 

00 10 K^l.IS 
IP-l 

Lll*IItK) 

L1=L1»L11 

L2=0 

DO 11 J>!ltIS 
NN^<K-1)*IS*J 
L22=I1 CJ) 

L2-L2+L22 

IRsO 

DO 12 JJ«ILiLl 

1R=1R*1 

1C*0 

DO 12 KK«IP»L2 
IC*1C*1 

12 A2tNN»IRiICI*AtJJiKK) 

11 IP«s1P*L22 

Ik«o 

DO 13 JJ«IL»L1 
IR«IR*1 

00 13 NK^liM 

13 B2JK,IH»KK)“a(JJ,KK) 

IL-lL+Lll 

10 CONTINUE 

RETURN 

END 



•• NO ERRORS** 
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COMMON ■ 02245 
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« 


(rZttO 

928 X 

0282 

DEM 

028* 

OEflS 

OZBfa 

0287 

0288 

0289 

0290 

0291 

0292 

0293 
029* 

0295 

0296 

0297 
0 Z 98 

0299 

0300 

0301 

0302 

0303 
030* 
03Qb 

0306 

0307 

0308 

0309 

0310 
03H 

0312 

0313 
031* 
031b 

0316 

0317 

0318 

0319 

0320 

0321 
D3Z2 
0323 
032* 
D3ES 


subroutine PPL(N,IW) 

c 

C POLE SHIFTING PROGRAM 
C 

’ integer do 

INTEGER YES 
REAL K1 

DIMENSION AiaO,10)fBLtlOt3)»A2(9,lDtIO)f82(3tlO,3)«nf5) •AtSiS), 
I01(b,S],QC3iS,S}f tlKSrS] ,B (5) ,K1 (6) ,R2{6] *0 (6) ,SK(6,10} tCOF (6) 

COMMON AlteiiAZ,B2,lI,A,0I,0,IU*B,RltR2>0,5K 
DATA YES/ZHYE/ 

URlT£nU,lS3) 

^53 FORWAT USX, "POLES OF THE SYSTEM") 

CALL ALAM(AfN,DfC0FiRltR2) 

DO 19 K=1»N 

19 MRITE(1M,100)K,RHXJ,H2(K) 

100 FORMAT (IXtflHLAHBDA *iFl*.6,*H * 

WRIIElIWfZOO) 

ZOO FORMAT llX ,"00 YOU WANT TO ALTER THE POLES 7, YES OR NOT") _ 
READ(1,330)OD 
330 FORMAT (AZ) 

IFtYES.NE.DDlGO TO 3*0 _ 

WRITEIIW»500) 

SOD FORMAT ax, "ENTER DESIRED EIGENVALUES") 

WRITE ClWtSlO) 

510 FORMAT ax, "ASSIGN INPUT UNITR/IX, "12") 

REAOtl,lll)lRO 
111 FORMAT (12) 

DO 209 J 1*1, N 
yRITEtlWiSOll 

SOI format llXf", _ t") 

209 REAOaRD«208}RUJl),RZ(JI) 

208 FORMAT (2F1*. 6) 

CALL DISP(H1,HZ,N,2) _ 

call KBAR(A,N,Z,IZ,D,B) 

DO 8 IslfN 

DO 8 J*1,N 

a AtI»J)*A(If J)*B(U*OtJ) 

WRlTEaWi9191 

919 FORMAT [ZOXfflHVECTOR K//) 

MHITEaw»I09) COlK),K*l,N) ' 

109 FORMAT (IX, 5F1*. 6) 

3*0 CONTINOb 
RETURN 
END 


Cl 

CJ 

<Y 

4 > 

« 

«3 

C> 
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PACE 0001 


FTN^^ compiler; HPZ4177 {SEPT-. 1974) 


0326 

0327 
Ct32t) 

0329 

0330 

0331 

0332 

0333 

0334 

0335 

0336 

0337 
0336 

0339 

0340 

0341 

0342 

0343 

0344 

0345 

0346 

0347 
0343 

0349 

0350 

0351 

0352 

0353 

0354 

0355 

0356 

0357 

0358 

0359 

0360 

0361 

0362 

0363 

0364 

0365 


SUBROUTINE AGRCIS) 

PfiOCHAM TO FORM AN AGGREGATE MATRIX 
REAL XB 

dimension A(10ylO)»BUOt3>tA2{9TlOflO)»B2(3ilDf3}rl! tS)TA4I5f5) 
1,A3(S»51 iQ 13, 5*5) till 15,5) ,B3 tS) ,RH6) »RZ(6) ,KBE6) ,SK(6tlO) , 
1C0F(6] *D(6),0l C5,5] 

COMMON A»a,A2,BZ,II,A4TA3,Q,Ill»B3,RltR2,K6,5K 

DO 5 K*ltIS 

Ll*lItK) 

DO 5 J«1,IS 
L2-IICJI 

IT*fK-l)*IS*J _ _ _ _ 

lFtK.EO.J)GO TO 13 
DO « L*1*L1 

DO 8 1*1,LZ _ _ _ 

8 A4Cl»U>A2(n»L»n 
DO 9 L-l.LZ 

DO 9 1*1, L2 . ... - . 

S>0. 

DO 10 IR*1,L1 

10 5»S*A4tlfIR)*A2aT»IRiL> , . 

9 aUl,L)=5 

CALL ALAH(Ol,L2tD*C0FpR2fKB) 

DO 11 I*1,L2 

IF(R2tn.LE.0>)I>0 TO 11 

R2(i)*saRTtH2an 

11 CONTINUE „ _ 

CALL BlOl (R2iLfL2) 

GO TO 14 

13 continue ^ - 

DO 20 1^1 fLl 

20 R2ll)»-5Xtl»Xl 

CALL SMAlHR2*L*L1) _ _ . 

R21U**K2CL) 

14 A3(K,J)*R2(L) 

5 continue . _ 

RETURN 

END 


** NO ERRORS** 
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PACE 0001 


FTN^^ COMPILERI HPZA177 tSEPT. 197A) 


0366 

1367 

0366 

0369 

0370 

0371 
03TH 
0373 
0376 

0375 

0376 
0J77 
03 7B 
03 79 
03B0 

0381 

0382 
03B3 
03U6 
03195 
0386 
0307 

0388 

0389 

0390 

039 1 

0392 
3393 
0396 

0395 

0396 

0397 

0398 

0399 
0600 
0601 
0602 
0603 
0606 

0605 

0606 
0607 
06.08 

0609 

0610 
0611 
0612 
0613 
0616 

0615 

0616 

0617 

0618 

0619 

0620 
0621 




SUBHOUTINE TRFIIS) 

PHQGRAM to transform SOaSYSTEM BY SIMILARITY TRANSFORMATION 


DIMENSION A(10ilO)iBaOf3] tA2(9t lOilO) iBZ C3«10t3) t II 15} f A1 C5t51 
li01(5t5}tqC3»5t5}tlll(5)fB3(5) 
lfC(2S)*LLt5)<HH[5) 

COMMON AtBtA2iB2in>AltQlf0iIllfB3 

00 25 K>lf£5 

LlsIItKJ 

00 26 3*1 tL I 

00 26 L>1*L1 


IF=U-lJ»Li*L 

26 cnP}«a(KfLtJ) 

call HiNVCCiLliDfLLtMM) 
DO 27 J=lfH 
DO 27 L»1*L1 
IP«U-1)»L1»L 

27 QULtJ]>CtIP) 

DO 33 JJ*lflS 

IT*JK-1)*IS*JJ 

L2*1IU3} 

DO 28 J»ltLl 

00 ZB L=liL2 
S>0. 

DO 29 IZ«l,Ll • 

29 5>S*aitJtlZ)«A2(lT»IZiU 

28 AltJ»U*5 
DO 31 J«liLl 
DO 31 L»lfL2 
S»0. 

DO 32 lZ*ltL2 

32 S6S*AltJfIZ}6Q(JJilZrU 
31 A2tITfJ»Ll=S 

33 continue 

DO 39 J-l.Ll 
S=0. 

DO 36 L*1*L1 
36 S>S*Qlt3iL}*S2tKiLfKi 
39 B3CJ1SS 

DD 20 J>ltLl 
20 B2(KtJiK}»B3{J} 

25 CONTINUE 
DO 72 M'ltlS 

DO 73 3»lfLl 
DO 73 L»1*L1 
lP*tJ-l)*Ll+L 

73 CClPJ=atiC»L»J} 

CALL KlNVCCtLliOfLLiHH) 
DO 74 

DO 74 L®1>U 
lP-tJ-lJ*Ll*L 

74 atKtLt3)>C(IPl 
72 CONTIN-Jt 

return 
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MCE BOOZ TRf FTN4 CDHPILERt HPZ4177 CSEPT, 197*) 
«*2Z END 
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0«Z3 

0%24 

D«»ZS 

Q^Ztk 

UAZ7 

OAZd 

0^29 

0A30 

0<»3l 

002 

0A33 

0434 

OOb 

043b 


SUBROUTINE ALAM(ATNfDiC0F*RltP2) 

C 

C PROCftAH TO calculate ElBENVALUES OF GENERAL N BY N MATRIX 
C 

OIHENSION ACSiS] »D(6) tCOF (6) «R1 (61 *R2(6) 

CALL COEFUAfNtDJ 

CALL POLRTlDiCOFiNiRltRZtlER) 

IFdERtEQ.OlGO TO IQ 
HKiTE(6ilOQ)lER 

loo FORMAT (IXtiOHERRQR CODE = *12) 

'10 CDNTINtlt 
RETURN 
END 


•» NO ERRORS** PROGRAM m Q006A 


COMMON * OOOOQ 


FACE 0001 


FTN4 COMPILERI WZM77 (SEPT. 


0M7 



SUBROUTINE UNITKRfN). 

04311 

C 



0439 

C 

SUBHOUTINE TO FORM UNIT MATRIX 

0440 

c 



0441 



dimension R(5tS) 

0442 



DO 9 K*liN 

0443 



DO 9 J>l.N 

0444 


9 

R(KtJ)a0. 

044<i 



DO a K«ltN 

0446 


S 

RCK*K)al. 

0447 



RETURN 

D44B 



END 

•* NO 1 

ERRORS** PROGRAM « 00072 COMMON a 00000 
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PACE 0001 


FTfiA COMPILER! HPZA177 {SEPT. 197AJ 


04A9 

SUBROUTINE COEFltAiN.O) 


O^SO c 



OASl C 

SUBROUTINE TO CALCULATE COEFFITIENTS OF CHARACTERISTIC 

POLYNOMIAL 

OAS2 C 

OF MATRIX A 


0<^S3 C 

’ 


OASA 

OIHENSION AlSib] tU t5i3} .Ct5>5} tO (6) 


O'tSS 

CALL UNIT1(B»MJ 


0<t‘i6 

DO 9 Kxl.N 


If 4b 7 

call PRODl(A.a»N»N.NJ 


OAbB 

CALL TRACI (B.NtS) 


04*39 

L*N-K*i , 


0460 

D(U»-I1./KJ*5 


0461 

CALL UNITKC.N) 


0462 

□1*0 il; 


0463 

call SCHHCfOl.N) 


0464 

9 CALL ADOKCfB.N) 


0463 

L»N*1 . 


0466 

DtL)«l. 


0467 

RETURN 


046S 

ENU 


«■ NO 

ERRORS** PROGRAM ■ OOEll COMMON * 00000 



I 

! « PAGE 0001 FTN^ COMPILER! HP2A177 (SEPT. 197A) 


0A69 SOeROUIlNE Sr.ALUtB.Ni FI 

OA70 C 

ovn C SUBROUTINE TO MULTIPLY MATRIX BY A SCALAR 
' O^TZ C 

0AT3 OIHENSION BtS) 

OAt*. ' DO a K»lfN 
UAT5 a B(X)>F«BtK) 

_ 0AT6 RETURN 

0A77 END 


•• NO ERRORS** PROGRAM » OOOAl COMMON a OOOOO 





PAGE a 001 


FTMA COMPILER* HPZ4177 CSEPT* 1974) 


SU4R0UT1NE VECPR(AfB)NiC) 

program to MULTIPLY MATRIX BY A VECTOR 

DIMENSION AtSiS)tBl5]*C{5) 

DO 9 K>IiN 
S*0. 

00 a J>1»N 
a 5>5*A(X«J)«BU) 

9 C(K)»S 
RETURN 
END 


A* NO ERRORS** PROGRAM ■ 0007S 


COMMON « 00000 


iO ^ 






ii 

li 

Ij 
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PACE SOOl 


FTN<k COHPILERJ HP24177, 


<SEPT. 1974) 


-y 

0517 

SVBROUrlNE PR001(Af8iN«ItH> 


05 la c 



0519 C 

SUaROOTINE 70 MUl.71PLy TNO MATRICES 

"S 

0520 C 



0521 

• OIMENSION AI5»5) *8 (5f 5) f C (5f5) 


0522 

S«0* ^ 


0523 

00 9 K>lfN f 


Q524 

DO 9 L«lfH 

- ■ — 

0525 

S«0. 


0526 

00 10 jslfl 

K* 

052^ 

10 S«S*A(KtJ) *81010 


0528 

9 C(K»L>*S 

• 

0529 

00 B KxltN 


0930 

00 8 L»1.H 


0531 

8 B(K(L>«C(KtU 


0SJ2 

RETURN 

- 

0533 

END 

■n 

•• NO 

ERRORS** PROGRAH > 00183 COHHON 


\ 


t 







I 


A.27 


PAGE 0<ei 


FTN* COMPILER* HP2A177 (SEPT. 197A) 


§534 


»MM0UTINE AE)OUA»BfN) 

§535 

C 

• 

0536 

c 

SUBROUTINE TO ADD TWO HA1 

0537 

c 


053§ 


dimension A(StS)*Bt5i5) 

0539 


DO S K^ltN 

- i §540 


DU • J«ltN 

054i 


B BtKf J)>b(KtJ)*A(KtJ) 

' “ 05*2 


RETURN 

- 0543 


END 


•• MO ERRORS** 


PROGRAM « 00066 


COMMON > OQOOO 





PAGE 9001 


FTN* COHPIUERI HPZhVn (SEPT. 1974) 





( 

i 




I 

t 

r 


PAGE QOOl FT^iA CCMPILERl MPZA177 «SEPT, i97*> 


0SS4 


subroutine SCHlCA(OtN) 

DSSS 

C 

• 

0556 

C 

SUBROUTINE TO MULTIPLY MATRIX BY A SCALAR 

0557 

c 


0550 


dimension A<5«5) 

0559 


00 9 Kal,N 

0500 


DO 9 J-IiN 

056 L 


9 A(KtJ)«U»AtKvJl 

056Z 


RETURN 

0563 


END 


•• HO «W<W5** PaOC’RAH ■ 00G57 COMMON « 00000 



* 






* 


# 


0 


0 


t 



1 1 




( 


H 

!i 



f 
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J 


1 


.1 


: ;l.v! j : . ■ J 


A.52 


A, 2, Optimization Program 

The only step that may introduce some ccmputational con^lexities in the 
optimization scheme described in Section 4 is the solution of the matrix Ric- 
cati equation for the evaluation of the local controls. Despite the fact that 
this computation is performed at the subsystem level and hence involves ma- 
trices of small orders, simulation on a digital con 5 >uter will invariably be 
necessary. Although many different methods for the solution of the Rlccati 
equation exist in the literature, the particular method that is adopted here 
is the iterative technique due to Kleinman [19] . In addition to determining 
the symmetric positive -definite solution P of the Riccati equation 

A^P + PA - PBR” Vp + Q * 0 , 

the program described here also conputes the eigenvalues of the matrices P 
and W = Q + PBR"‘‘‘B^P that is necessary in the evaluat5.on of the suboptimal- 
ity index e . 

The simulation analysis was conducted on the HP 2100 digital conputing 
system (32K memory) in FORTRAN language. In the following description, only 
the subroutines MINV, SMJ and POLRT are to be sipplied externally (from IBM 
Scientific Subroutine Package), ’idiile the rest are contained internally. 

Since the computation involves only the subsystems that result from a su5.ta- 
ble decomposition of the overall system and hence are necessarily of small 
dimensions, the program is prepared to handle subsystems of dimension up to 
five, 

PffSfffflPTIOff OF THE EXTERML SVBHOUTimS (From SSF) 

Subroutine MINV 
PMrpoae; 


Invert a matrix. 


JJsage: 

CALL MINV CA, N, D, L, M) 

Oescvi'ption of parameievs: 

A - Input matrix, destroyed in computation and replaced by 
resultant inverse. , 

N - Order of matrix A . 

D - Resultant determinant, 

L - Work vector of length N . 

M - V/ork vector of length N , 

Remxpks: 

Matrix A must be a general (nonsingidar) matrix, 

Subroutims and function sidiprofframs required: None. 

Method: 

The standard Gauss -Jordan method is used. The determinant 
is also calculated. A determinant with absolute value less 
than 10 (-20) indicates singularity. 

Subroutine SIMQ 
Purpose: 

Obtain solution of a set of simultaneous linear equations AX ~ b 
Usage: 

CALL SBJQ (A, B, N, KS) 

Description of parameters: 

A - Matrix of coefficients stored solumnwise. These are des- 
troyed in the computation. The size of matrix A is N 
by N . 

B - Vector of original constants N) . These are re- 

placed by final solution values, vector X . 


A. 54 


N - Number o£ equations and variables. N must be greater 
than 1. 

KS - Output digit: 0 for a normal solution j 1 for a sing- 

ular set of equations. 

Remccpksi 

Matrix A must be general. If matrix is singular, solution 
values are meaningless. 

Subroutines and function subprograms required: None, 

Method: 

Method of solution is by elimination using largest pivotal divisor. 

Su broutine POLRT 
Purpose: 

Confutes the real and coniplex roots of a real polynomial. 

Usage: 

CALL POLRT (XCOF, COF, M, ROOTR, ROOTI, lER} 

Description of parameters: 

. XCOF - Vector of M+1 coefficients of the polynomial ordered 
& 

from smallest to largest power, 

COF - Working vector of length M+1 . 

M - Order of polynomial. 

ROOTR - Resultant vector of length M containing real roots of 
the polyncmual. 

ROOTI - Resultant vector of length M containing the correspond- 
ing imaginary roots of the polynomial. 
lER - Error code where 

lER * 0 No error 

lER * 1 M less than one 


lER = 2 M greater than 36 
lER = 3 Unable to determine root inth 500 
iterations on 5 starting values, 
lER = 4 High order coefficient is zero, 

Remaz^ks: 

Limited to 36-th order polynomial or less. Floating point over- 
flow may occur for high order polynomials but will not affect the 
accuracy of the results. 

Subroutines and function aidtprogrcans required! None 
Method: 

Newton-Raphson iterative technique. 


PACE OUDl 


I 


COHPILERt HP2417" JSEPT* 1974) 



UUCl 

0062 PtiOGSA*! «1CAT 

OoOi C. SOUOTlOH OF ALGEBHftlC RICCAJI EOUaT 10» BY rtLEINMAN ITERATION 
OIHEMSIUM AI£i.b)tB(S«5))atSf5)fR(St5}tSU(5iB)TS(5t5)t 
OUCb iBC(aiS) iAM(b*b) jTStifS) fAUXtSfS) .P'S*5» *Wl5»£iJ t 

OOUb 2KIMC»51*5F lS.3)tUVlb).HVi:S)*VRl2B)»TB(Si5)»TA{St5)» 

DOOT CAilAilBtb) lAuxbtBiB] tAOXClStSl (AUADtCfS) tAUAE[St5) tAUXFtB.S) > 

OOtrd ‘»auA(j|StB}iAUXH(btb)«VlM6> iCOF 16) (6J t«S16) 

OoO-^ C htAU Tit STATE OlFithSlON N 

OOlU Rt.Hj|b«nN 

UOll 1 FUSMATUb) 

ttH<r L F-r.ilU TMt CUNTROL DIPt^5!0^l H _ . 

OOlJ HtflUtbllJH 

ocm C RtAO IHL SYSTEM HAIKU A 

OJlb hEAUl6*2j C{A(lfJ)»J=X*NJ»l=l.N) 

Quid 2 ) UKHAT(t)Fl(i«b) 

OOH C KtAU the STmTE-CUNTHOL HATRIA B 

oOla ttEAOlb.2) ItHII.Jj»J=liH)tI=lfN) . . . „ . _ 

OSJlY C KtAu f*(t CObl HAtKICEb u AhD H 

uU2U Ktni)(bt’2] J}fJ=ltN) fUlfN) 

0021 RtAiJ(&,2) n«Cl»JJ *I=XrH) ... , . . 

0022 KHlIE(bt20J 

U02J 20 FORMAT (JK*"b(lLim UN OF KlCCATl EQUATION-KUEINHAN ITRATIVE METHOD’* 

002** . . 

'JU2s WMi1EI6»2i) 

002b 21 FortMAT(bX»"SYbTtM HAlHlX A"t//) 

0027 DU 22 I=1.N .. 

U02d 22 Y>KI(£(bt23) tACliJ) 

UU2Y 2J FORMAT UHtitlU(2A«E13.S] ) 

0040 *.H|TE(b»24) _ 

Ouil 24 FUR*iat(//2) 

0«J2 WRITE (b«2b) 

OOiJ 2b FORMAHbXi'iSTATfc“COHmOU MATRIX 

OOJ** DO 2b 1 = 1 »h 

tiojb 2b wRlTEtbt233 (B(ltK) tK=liM) 

yo4b wKirtt6»241 ,, _ _ 

DoJ7 w**|(t(6)27) 

004b 2T FOHMATIbA.t'CUbl MATRICES 0 AND 

0047 wRilECb>2B) . . . 

0040 2B F0RmATURQ»"0 HATHIX'i«/2) 

UU**1 UU 2T 1 = 1 fN 

o04^ 29 bKnE<bi23Mun(J]tJ=ltN) 

0043 icill£(b>30) 

004H 40 FOmMAT U/ ilXi"R MATKlX“i//) 

004b UU Jl l=liH . . 

0U4b 41 hkITF tbi23) CKIli Jl iJ=liM) 

0C47 KK|;t(bl24) 

Oo4tt c RtAU TnE 41ARTIN0 APPKOAlMANT SO 

0W47 wLAlMbiO ICbU(liJ>.J=liN}fl=l,H) 

OubO no 3 1=1 tM 

OObl DO 3 J=ltN 

0Db2 4 btl ij)ebt'tl« J) 

U0b4 CAU. «*uvitHiM,VK) 

rf)b4 C“',L Ni\\/{vR»H*utt.LV»'4V| 

OObs CAU. YlP4IV«tM*HlN) 

UObb t FHiT 1 ttN CUMMENCtS here 



r\ 


PAGE 0D02 HICAT FTN* COMPILER: HP24177 (SEPT. 197«) 


kRlTE{6f33) 

FOKi^AT tlOAt«STAHT OF THE FIRST ITERATION''.///) 
ii-<iTE(5,3'*) 

FORMAT (bX."STARTlNG APPRDXIHANT-S MATRIX".//) 

00 J5 1=1»M 

Hrfn£It>t23) tS(I.J) . 

•fril Tt. (6.E**) 

DO £00 nVAH=l.£>0 

call >:AMUL(B.S»b5»NtM.fl3 . — . ... .. 

OU 4 I=1.N 
00 4 J-liN 

0Sa.J)=-b3(I»J) . - 

THE M03IF1EU srsTfcH MATRIX AH IS OBTAINED NOW 
call «»*SUH(A.kJ5.AH.NfN) 

WR1[E[6.36) . . 

FORMAT tbX."MUDlFI£D SYSTEM MATRIX— AMrA-B*S"t//) 

BO JT 1*UN 

RrflTE(iit23) (AM(1»J).J*1.N) 

riRlTE((>.2A) 

CALL ElVAL(rtH,N.VD.C0F.Hl.R2) 

KrclfLtbtllB) 

FORMAT tbX."UOENVALOES DF CLOSED LOOP MATRIX".//) 

00 1A3 J-ltN 

WMlIttf.t2i)RUJ).KElJ) . . 

THE L7APjf,tn# MATRIX EQUATION IS NOW FORMED 

call MBTRh(SiTS.MtN) 

call MAHJL(TS»R»AUXb.N.M.H) 

CALL HAHUL (AUXbtSiAUXHiNtH.N) 

call MAaUM(U.AUxH.D»M.N) 

TnE LYAPOMUV HAIHIX EQdATIDN IS NOW P«AH*TR(AHJ "P»- 
CALL LyAPUtU.AM.P.N) 

WMI1£(5)B) 

FORMAT ISX."P HATKIx-SOLOTION OF LYAPUNOV EQUATION". 
DO V Isl.N 

■rilTEtfiflO) IPtt.Jl .J-ItNl . „ 

FORMAT (lHa.lO(2XfElJ*S) I 
KRl TE(&>2A) 

compute the APPrtOXlMA>NT FDR THE NEXT ITERATION 
CALL MATKNtri.TUiN.H) 

CALL MAHOLCRIN.IU.AUX.H.H.N) 

CALL HAvuLCAOX.P.SF.M.N.N) . „ _ 

CALL MAjRNtA.TA.N.Nl 
Call MAHULtTA.P.AUXA.N.NfN) 

CALL MAMJLtP.A.AUXIj.N.N.N) - 

CALL HAMULlb.SF.AUXC.N.M.N) 

CALL MAMliLtPiAUXC.AUX.N.NcNJ 

COMPUTE The a MATRIX . W-Q*F»B*RIN*T3*P 

CALL HASUH(U.AOXtWtNtN) 
a-miTE Ct.l31» 

F'UrtHAT CSX»"« MATRIX WAB.P*B«RIN*‘Ta*P".//T 

UU 132 1*1 fN 

’ •.RnEt&»23nHtI.^}.0=liN} 

VRltEtS.OA} 

DO 112 1=1. N * 

GO 112 J=l.w 
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HAUE 


G003 RlCAT FTN-i COMPILER* MP24177 (SEPT. 1974) 


0LI3 

dll'* 


112 

Ollb 

0116 



OUT 

Olio 

C 


0119 


111 

0120 


lU 

0121 

0122 


ll4 

0123 

Dl2" 


11b 

Ui2t> 

c 


012U 

c 


I'lZT 

0120 

0)29 

. c 

- 

DUO 


0 

OUl 

c 


0132 

c 


UUJ 

DU*. 

Ul3b 

c 


0136 

013/ 


133 

ulJd 

0139 

omv 

0191 


134 

VI92 

Ol4j 

-13b 

OImh 


136 

0146 

0196 



014/ 

OlAB 


3b 

Ol49 

UlbO 


106 

Olbl 

Dlb2 


39 

OlbJ 


200 

CrtS4 

Olbu 


99 


AUXtI»J)»-AUMI.J) 

CALL HASOHIAUX.AUXAiAUAF.NiNl , 

CALL HAStJM(AUXF»AUXLtiAiJXD»NtH) 

CALL HASUM(AUAU»U>AIiaE*NiN) 

ChLCKlNb THE SOLUTION 

WKir£Js,un , 

FOKrtAT(30Xt"CHELKlNlj THE SOLUTION”*//) 

WHITE (6*11*.) 

FOH«AT tbA*”VALUc. OF KlCCATl MATRIX".//) . ^ . „ 

00 115 1=1. N 

WRITE t&*23] (AUXE(l.J) . JsltN) 

WHlTEtb.SA) _ 

matrix SF is the SIAKTINC APPROXIMANT FOR NEXT ITERATION 
ThE STEADY VALUE OF SF IS THE FEEDBACK GAIN MATRIX FOR THE OPTIMAL 
KtbULArufi K=K1N*1U«P 
00 9 I=1.K 
DU b J=1*N 
stitj)=bFn.j> 

AUUITIONAL commutations— CALCULATION OF EIGENVALUES OF 
HAThICES W and P— This PART MAT BE OMITTED WHEN THE SOLUTION 
Of A HlLCATl EUUATION ONLY IS DESIRED . . 

CALL EIVAL(P.N*VD.C0F»Rl»R2l 
WH1T£(5.133) 

FORMAT tbX."ElG£NVAUUES OF P”) i//J . _ „ 

DO ilA Jsl.N 
WRITE t&.EJlRl CJ) *R2( J) 

WMilEtO.EA) . . 

ElVAL<W*NtVU.C0FtRl)K2) 

WKir£(6.135} 

FORMA! tux. "EIGENVALUES OF W"f//I _ „ 

no 116 1 = 1 iN 

WRITEt6t23)Rl(l).R2U) 

wRirEte.EA) .. _ 

WRITE [6.3HI 

FUhMATtlOX."NEAT ITERATION COMMENCES HEBE".///) 

»hXIE(6.1Q6) . . , _ , 

FORMAT (bx , "modified APPROX INANT-NEW S MATRIX".//) 

DO 39 1=1. H 

WRITE 16.23) (5U*J).J»1.N1 . _ 

WRl1E(6.24) 

CONTINUE 

STOP . _ _ 

END 


NO ERRORS** PROGRAM * 02921 


COMMON s 00000 


^ fA&e UOOl FTN*f -COHPILERI HP24177 (SEPTi 197A) 



Ulbb 


SUUKOJTINE HAM0L(F»G.H.N1,H2,N3) 

4 

0167 


OIHENSIUN F(S,6).G(6.S)iHC5i5) 


Cibd 


00 51 I=1.N1 


. 0169 


00 61 J=1.N3 

\ 

UlbU 


SUM^O.D 


Olbl _ 


DO 62 i(=l]N2 .. 

1 

^ OloZ 

52 

SU»i=SlM*r ( 1 .KJ »G (N. Jl 

i 

01b j 


imij]=bUM 

\ 

' Oltj** .. 

. 51 

CUNllNOt - ^ . 


^ 016b 


RETURN 

t 

OtbQ 


END 


•* NO E(WORS»* PflOOHAH = OOD91 COKHDN * 00000 


t 

< PAGE 0001 FTNA COMPILER! HPZU77 (SEPT. 197«) 

• 

) 

, " Oib/ ' SUdRUUUNE HASOM(f.GTHtNltN2} ’ ” 

‘ - OlbB DIHLNSION F tS.bl f&t5»S> fM(Eif5) 

’ „ 01t>9 _ OV 01 I=1»NX . _ 

’ sure 1)0 61 j=i«d2 

• - Oin HdtJlsF 

_ tfltd 61 COiTlKUt _ 

01/J reiukn 

01/** ENG 


*• NO CRH0HS«* program = 00067 



COMMON ■ 00000 


PAGE 0001 FTNA COMPlLEHt HP2A17T (SEPT. 197^1 


Dl/b SUaKQUtlNE MArRN(FiGfNl*M21 

Vl/G UlP-tNSIDN F(Stbl iGtbtS} 

Ul// OU 71 1=1. N1 

017a Ou 71 J*1.N2 

Cl/9 G(J.T)zF(ltJ) 

Cl6t) . 71 CUNtIWiJt. 

Olal KETURN 

GIB2 E«G 


*• NO errors** program = OOOE5 COKhON * OOOQO 
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.w 
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SUSKOUTINE LYAPUCG»A.riiM) 

DIHtNSIUi^ G(5t5)iA(bfS} »rtl5»5>»Pa5) iLCStSl *Ua5*15} i 

M=tti*i}*N/2 

CALL T'?iM(G»P,N»rt} 

CAUL MATHHtMtLJ 

. CALL LP(L«A»NiMtU> , 

00 in I=1»H 

10 PU)!=-PCI) 

CALL HAVECtU»H*PU _ - 

CALL SlMaiPltPiHtKSJ 

. 12 HRlIEibirlOO) . 

100 F0kMATMA»13HSIN£iULAtt CASE) 

11 CAUL ATHI(PtH*H»NI 

retohn. . 

EMU 


NO EHROkS** 


PROGRAM 


COMMON X 00 000 


PAGE 0001 


FIN4 compiler* HP24177 (SEPT. 1974) 


SUtJKOUUNE TklN(Aip.NiM) 
QIHENSIUN ALS.bliPaS) 
M=(l . .... 

DU 6 Kxl.N 
L=K 

DU-d X»LiH 

H=M*1 

fc P(M)sA(L.n 
return - 

END 


** NO ERRORS** 


PROGRAM X D0064 


COMMON X 00000 
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PACE OltOl 


FTN‘^ COHPILERI HP24177 tSEPT. 1974) 


0<dOV S’JBKPUTINE ATPI{P.AfMf"l) 

J<!n) OIMEMSION ACbtb)»P(i;>) 

Q<ill MaU , 

uaia DO a K=itN 

oai3 

. _ Uai**. __ DO t# l=UiN 

UaiS K=M*1 

oatb a A(U«I)aP(M) 

0ai7 . , DU V J=IiN .. ... 

caia L=j 

Oai9 00 9 I=L»N 

adat) . 9 AUiji=AijfU „ 

tidal RtTUHtj 

odaa . Enj 


'•* NO ERKOHS** PROGRAM = Q0l04 COMMON « OOOOO 

PAbE OUOl FTN4 COMPILER* HP24177 <SEPT, 1974) 



0223 SOaROuriNE MATLCNaKiL) 

0224 DirttNSION H5t5)tLl(lS) 

oaap . . DO B K=i,H 

022b a Ll*X)=K 

022* MlaO 

. . U22B DO IB <=I»N «. _ 

0229 , Lb=K 

D 2 iO 00 la l=LSfH 

D2J1 

0232 Vi LlLb,i)aLHMl) 

Q2J3 00 9 J=ltN 

. - 0234 LS=J 

023b DO 9 I=L5iN 

0233 9 LUtJ)-LtJ»l) 

023/ RETOKti 

OdJa ErtU 


'** hO ERftORS**" PROOHAM « 00128 .COMMON « OOOOQ 

page 0001 FTN4 COMPILER* HP24177 .(SEPT. 1974) 


i 

0239 

SUBRUUnNE kAvecu«h<p)} 

— ■ 

oa^o 

UlMtNSJON AdSiIb) «P< (sas) 


. 02I* L 

DU t K-l^K 


,024a 

1)0 a J=liN 

'■ 

0243 

L=N*tK'l)*3 


0244 

a pilL)=AtJ*h) 


024a 

RETURN 


U246 

End 


NO ERRORS* 


"ROORAM * 00061 


COMMON s OOOOO 
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» 
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SUoHOUTlME MAV£(A.M.P1J 
UIMfcNSIONI A(5»S>»PH25) 
DO a.K=l,N 
DU a JsltS 

(S PllLl*AtJfKl 

RETUBiK 

EMU 


•* NO ERRORS** PROGRAH « 00061 


COMMON > 00000 


PAGE 0001 


FTN*. compiler: HP24177 (SEPT, 1974} 


0H55 SUtsHOUTlNE V£HAtPl,NiAJ 

02bo DIMENSION AtStbJ ,P1 (251 

0:ib/ . „ ...Do tt K=ltN _ . 

Odad 00 & JxltN 

02bV L«N*tK-l)*J 

0260 B A(J»RJ=P1(L) . . . 

0261 RETORN 

0262 END 


♦* NO ERRORS** PROGRAM » 00062 


COMMON ■ OOOQO 


PAGE 0001 


FTN4 compiler: HP24177 (SEPT. 1974) 



D26J 


SubROUTlNE U 


0264 


DlP-LNSIUN L(; 


. 0£6b 


OU 11 1=1 iM 

• 

' 026b 


DO 11 J=i»M 

^ ■* 

0261 

11 

uap3}=u. 


02&6 


DO 12 1=1 »N 


0269 


DU 12 J=1*N 


02/0 


Oo 12 <=1,N 


ocn , 


Il=Ltli2J 


■ 0272 


IJ=LIJ,M 


0273 

12 

U(Il.lJ]i:A(Ji 


.. 0274 . _ 

1. *. .»• 

DO 13 1=1, D 

' 

I>2/6 


DO 13 J=1»M 


0276 


iL=LiiTn 


0277 

13 

OtlL.3)=2.*U 

' 

0276 


RETURN 

Mi*- 

0279 


ENU 

m' 

•P NU 

ERHO.iS** PSDGKi 


COMMON s ODCiOO 


OEIGIWMI tAGE IS 


q<:bi 

JldUi 

Q<'B4 

Vii6b 

vdar 

v^tia 

OifO^ 

o^yo 

c 


** Nt» 


oa?i 

tiiiva 

bayj 

02^4 • 

Vd’tb 

fli’ff, „ 

Qd'tT 

oavd 

oa9» 


.*■•. -.NO. 


0300 

0301 

0302 .- 
uJ03 
03D4 
0305 _ 
CtiOb 
030/ 

f » 

'.t 

•’-_0 

(Jill 

0312 

OJlJ 

0314 

0Jl5 


** NO 
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SOBKOUTINE EIVALtAiNiOtCOF.Hl *«2) ‘t 

OlNtNSION 4tSib},D(6J,CC!Ft6),Rl(6)fRe(6) ‘ 

... CALU COfcFUAtNiD) .-. _ - 

Hi*N*l 

CALL P0URTttltC0F»N*Rl»R2»ieRJ 

IFtlER.EQ.OOO TO IQ „ .... ' 

WrtlTElb.lOOUEH ft 

lUQ FOrtHATtlXilJHtHKOH C00£ * .12) 

„.xa CONTtKUfc. . . .... - . * 

RETORN . 

END 


EttkORS** PROGRAM = Q0064 COMMON = OODOQ 
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SUdNOUTlNE UNlTHRiN) 
UXMbNSION RtSiS) 

- UO 9 N=1.N 

00 9 J=liN 
9 H(K*3)«U. 

. DO a K*liN 

tt R(K*K}«1. 

RETuHN 

ENU 


ERRORS** RRD6RAtt..=_(U3072 COMMON =_0 0000 
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SOUhoUTINE COEFl «A,W,U) 
dimension AtStS) tblS.S) tCOtS) |D{6) 

CALL UNITKOiN) .. „ 

DO V K=1.N 

CALL PRODlCAiB.NtNiN) 

C«LL TRACI tb^NcSl .. . 

L=N-.tf 1 

D(L>=-a./K)*S 

_ . CALL ONITl (L.N) . . . . 

01=0 (U 

CALL SCiiKCiOliN) 

.. 9 CALL A3OUCj0iN) . _ 

L=H*1 

OtD-l. 

HfclOSK _ 

ENt) 


ERRORS** program a 00211 COMMON *00000 



a I ^ 

• . 

• 4s. 

W 
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SubKOUTl^E liuC«{A.DtN} 
DlMEHSlUN A 5tbJ 
Oti ^ S=1.N . 

(lU 9 J=1«N 
9 J)eu*A(K.J) 

RLTUKM 

Er^D 


MU ERKOPS’** PROGRAM = OOOS7 COHMOK = OQflDQ 
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SU3B0UTINE THACHA.N.Sl 
DlHtNSIOfS At5»5) 

S=tl. 

UO ti K=1.M 
B S=S*A(AtK) 

RtTUMN 

ENU 


NO ERRORS** PKOGHAH = DOOAJi 
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SUBROUTINE ADOKA.bfNt 
DIMENSION AvSiS) «Ut5«5) 
DO a K=lfN 
DU o JxltN 

a atKtJ)sa(KTJ}«A(K*J) 

RtIURN 

FNU 


** NO ERRORS** 
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PROGRAM s 00066 


COMMON « 00000 


FTnh^ compilers HPZ4177 (SEPT. 197AJ 


SUOKOUnNE PROJKA.atNtItM} 
DIMENSION AtS«b}tBlSt5|.CCb*5) 

s=o . 

, DU 9 R=1,N 
DU 9 L=i»K 
S=0. 

UU IQ J=Itl 
ID S=S*A(H,JI«GtJ,L» 

V C(AtLl=S 
DO a K=IiH 
DO a L-l.H 

a aiA»L)=c(N«u . _ 

KLrUKN 

END 


NO tPROrtS** PROGRAM = 00103 COMMON * 00000 
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bJbl End A 


